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We describe azimuth structure commonly associated with elliptic and directed flow in the context 
of 2D angular autocorrelations for the purpose of precise separation of so-called nonflow (mainly 
minijets) from flow. We extend the Fourier-transform description of azimuth structure to include 
power spectra and autocorrelations related by the Wiener-Khintchine theorem. We analyze several 
examples of conventional flow analysis in that context and question the relevance of reaction plane 
estimation to flow analysis. We introduce the 2D angular autocorrelation with examples from 
data analysis and describe a simulation exercise which demonstrates precise separation of flow and 
nonflow using the 2D autocorrelation method. We show that an alternative correlation measure 
based on Pearson's normalized covariance provides a more intuitive measure of azimuth structure. 

PACS numbers: 13.66.Bc, 13.87.-a, 13.87.Fh, 12.38.Qk, 25.40.Ep, 25.75.-q, 25.75.Gz 



I. INTRODUCTION 

A major goal of the RHIC is production of color- 
deconfined or QCD matter in heavy ion (HI) collisions, a 
bulk QCD medium extending over a nontrivial space- 
time volume which is in some sense thermalized and 
whose dynamics are dominated in some sense by quarks 
and gluons as the dominant degrees of freedom [l| . "Mat- 
ter" in this context means an aggregate of constituents 
in an equilibrium state, at least locally in space-time, 
such that thermodynamic state variables provide a nearly 
complete description of the system. Demonstration of 
thermalization is seen by many as a necessary part of the 
observation of QCD matter. 

A. Global variables 

One method proposed to demonstrate the existence of 
QCD matter is to measure trends of global event vari- 
ables, statistical measures formulated by analogy with 
macroscopic thermodynamic quantities and based on in- 
tegrals of particle yields over kinematically accessible mo- 
mentum space. E.g., temperature analogs include spec- 
trum inverse slope parameter T and ensemble-mean p t 
f>t- Chemical analogs include particle yields and their 
ratios, such as the ensemble- mean K/tt ratio Q. Corre- 
sponding fluctuation measures have been formulated for 
the event- wise mean p t (pt) ("temperature" fluctuations) 
and K/tt ratio (chemical or flavor fluctuations) [1, 0, @ • 
Arguments by analogy are less appropriate when deal- 
ing with small systems ('small' in particle number, space 
and/or time) where large deviations from macroscopic 
thermodynamics may be encountered. 

B. Flow analysis 

One such global feature is the large-scale angular struc- 
ture of the event-wise particle distribution. The compo- 
nents of angular structure described by low-order spheri- 



cal or cylindrical harmonics are conventionally described 
as "flows." The basic assumption is that such structure 
represents collective motion of a thermalized medium, 
and hydrodynamics is therefore an appropriate descrip- 
tion. Observation of larger flow amplitudes is therefore 
interpreted by many to provide direct evidence^ for event- 
wise thermalization in heavy ion collisions [f|. Given 
those assumptions each collision event is treated sepa- 
rately. Event-wise angular distributions are fitted with 
model functions associated with collective dynamics. The 
model parameters are interpreted physically in a thermo- 
dynamic (i.e., collective, thermalized) context. 

However, collective flow in many-body physics is a 
complex topic with longstanding open issues. There is 
conflict in the description of nuclear collisions between 
continuum hydrodynamics and discrete multiparticle sys- 
tems which echoes the state of physics prior to the study 
of Brownian motion by Einstein and Perrin 0, |8] . Be- 
yond the classical dichotomy between discrete and con- 
tinuous dynamics there is the still-uncertain contribu- 
tion of quantum mechanics to the early stages of nuclear 
collisions. Quantum transitions may play a major role 
in phenomena perceived to be "collective." Premature 
imposition of hydrodynamic (hydro) models on collision 
data may hinder full understanding. 



C. Nonflow and multiplicity distortions 

A major concern for conventional flow analysis is the 
presence of "nonflow," non-sinusoidal contributions to 
azimuth structure often comparable in amplitude to sinu- 
soid amplitudes (flows). Nonflow is treated as a system- 
atic error in flow analysis, reduced to varying degrees by 
analysis strategies. Another significant systematic issue 
is "multiplicity distortions" associated with small event 
multiplicities, also minimized to some degree by analysis 
strategies. Despite corrections nonflow and small multi- 
plicities remain a major limitation to conventional flow 
measurements. For those reasons flow measurements in 
peripheral heavy ion collisions are typically omitted. The 



2 



opportunity is then lost to connect the pQCD physics of 
elementary collisions to nonperturbative, possibly collec- 
tive dynamics in heavy ion collisions. 

D. Minijets 

A series of recent experiments has demonstrated that 
the nonsinusoidal components of angular correlations at 
full RHIC energy are dominated by fragments from low- 
Q 2 partons or minijets [1 0, El El ■ Minijets in RHIC 
p-p and A-A collisions have been studied extensi vely via 
fluctuations [j| [l(| and two-particle correlations [ill . 
Minijets may dominate the production mechanism for the 
QCD medium [l3|, and may also provide the best probe 
of medium properties, including the extent of thermal- 
ization. Comparison of minijets in elementary and heavy 
ion collisions may relate medium properties and collec- 
tive motion to a theoretical QCD context. 

Ironically, demonstrating the existence and properties 
of collective flows and of jets (collective hadron motion 
from parton collisions and fragmentation) is formally 
equivalent. Identifying a partonic "reaction plane" and 
determining a nucleus- nucleus reaction plane require sim- 
ilar techniques. For example, sphericity has been used to 
obtain evidence for deformation of partic le hp t/Et angu- 
lar distributions due to parton collisions [l4( and collec- 
tive nucleon flow [l^]. At RHIC we should ask whether 
final-state angular correlations depend on the geometry 
of parton collisions (minijets), N-N collisions or nucleus- 
nucleus collisions (flows), or all three. The analogy is 
important because to sustain a flow interpretation one 
has to prove that there is a difference: e.g., to what ex- 
tent do parton collisions contribute to flow correlations 
or mimic them? The phenomena coexist on a continuum. 

To resolve such ambiguities we require analysis meth- 
ods which treat flow and minijets on an equal footing 
and facilitate their comparison, methods which do not 
impose the hypothesis to be tested on the measurement 
scheme. To that end we should: 1) develop a consistent 
set of neutral symbols; 2) manipulate random variables 
with minimal approximations; 3) introduce proper statis- 
tical references so that nonstatistical correlations of any 
origin can be isolated unambiguously; 4) treat azimuth 
structure ah initio in a model-independent manner using 
standard mathematical methods (e.g., standard Fourier 
analysis); and 5) include what is known about minijets 
( "nonflow" ) and "flow" in a more general analysis based 
on two-particle correlations. 



E. Structure of this paper 

An underlying theme of this paper is the formal re- 
lation between event-wise azimuth structure in nuclear 
collisions and Brownian motion, and how that relation 
can inform our study of heavy ion collisions. We begin 
with a review of Fourier transform theory and the rela- 



tion between power spectra and autocorrelations. That 
material forms a basis for analysis of sinusoidal compo- 
nents of angular correlations in nuclear collisions which 
is well-established in standard mathematics. 

We then review the conventional methods of flow anal- 
ysis from Bevalac to RHIC. Five papers are discussed in 
the context of Fourier transforms, power spectra and au- 
tocorrelations. To facilitate a more general description of 
angular asymmetries we set aside flow terminology (ex- 
cept as required to make connections with the existing lit- 
erature) and move to a model-independent description in 
terms of spherical and cylindrical multipole moments. We 
emphasize the relation of "flows" to multipole moments 
as model-independent correlation measures. Physical in- 
terpretation of multipole moments is an open question. 

We then consider whether event-wise estimation of the 
reaction plane is necessary for "flow" studies. The con- 
ventional method of flow analysis is based on such esti- 
mation, assuming that event-wise statistics are required 
to demonstrate collectivity, and hence thermalization, in 
heavy ion collisions. 

We define the 2D joint angular autocorrelation and de- 
scribe its properties. The autocorrelation is fundamental 
to time-series analysis, the Brownian motion problem and 
its generalizations and astrophysics, among many other 
fields. It is shown to be a powerful tool for separating 
"flow" from "nonflow." The autocorrelation eliminates 
biases in conventional flow analysis stemming from finite 
multiplicities, and makes possible bias-free study of cen- 
trality variations in A-A collisions down to N-N collisions. 

"Nonflow" is dominated by minijets (minimum-bias 
parton fragments, mainly from low-Q 2 partons) which 
can be regarded as Brownian probe particles for the QCD 
medium, offering the possibility to explore small-scale 
medium properties. Minijet systematics provide strong 
constraints on "nonflow" in the conventional flow con- 
text. Interaction of minijets with the medium, and par- 
ticularly its collective motion, is the subject of paper II 
of this series. 

Finally, we consider examples from RHIC data of au- 
tocorrelation structure. We show the relation between 
"flow" and minijets, how conventional flow analysis is 
biased by the presence of minijets, and how the autocor- 
relation method eliminates that bias and insures accurate 
separation of different collision dynamics. 

We include several appendices. In App. A we review 
Brownian motion and its formal connection to azimuth 
correlations in nuclear collisions. In App. B we review 
the algebra of random variables in relation to conven- 
tional flow analysis techniques. We make no approxima- 
tions in statistical analysis and invoke proper correlation 
references to obtain a minimally-biased, self-consistent 
analysis system in which flow and nonflow are precisely 
distinguished. In App. C we review the mathematics 
of spherical and cylindrical multipoles and sphericity. In 
App. D we review subevents, scalar products and event- 
plane resolution. In App. E we summarize some A-A 
centrality issues related to azimuth multipoles and mini- 
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jets. 



II. FOURIER ANALYSIS 

The azimuth structure of nuclear collisions is part of 
a larger problem: angular correlations of number, p t and 
Et on angular subspace (771,772,01,02)- There is a for- 
mal similarity between event-wise particle distributions 
on angle and the time series of displacements of a parti- 
cle in Brownian motion. In either case the distribution 
is discrete, combining a large random component with 
the possibility of a smaller deterministic component. The 
mathematical description of Brownian motion includes as 
a key element the autocorrelation density, related to the 
Fourier power spectrum through the Wiener-Khintchine 
theorem (c/. App. |A"|). 

The Fourier series describes arbitrary distributions on 
bounded angular interval 2ti or distributions periodic on 
an unbounded interval. The azimuth particle distribu- 
tion from a nuclear collision is drawn from (samples) 
a combination of sinusoids nearly invariant on rapidity 
near midrapidity, conventionally described as "flows," 
and other azimuth structure localized on rapidity and 
conventionally described as "nonflow." The two contri- 
butions are typically comparable in amplitude. 

We first consider the mathematics of the Fourier trans- 
form and power spectrum and their role in conventional 
flow analysis [lg|. We assume for simplicity that the 
only angular structure in the data is represented by a few 
lowest-order Fourier terms. In conventional flow analysis 
the azimuth distribution is described solely by a Fourier 
series, and corrections are applied in an attempt to com- 
pensate for "nonflow" as a systematic error. We later re- 
turn to the more general angular correlation problem and 
consider non-sinusoidal (nonflow) structure described by 
non-Fourier model functions in the larger context of 2D 
(joint) angular autocorrelations. Precise description of 
the composite structure requires a hybrid mathematical 
model. 

Event-wise random variables are denoted by a tilde. 
Variables without tildes are ensemble averages, indicated 
in some cases explicitly by overlines. Event-wise averages 
are indicated by angle brackets. The algebra of random 
variables is discussed in App. [B] Where possible we em- 
ploy notation consistent with conventional flow analysis. 



A. Azimuth densities 

The event-wise azimuth density (particle, pt or E t ) is 
a set of n samples from a parent density integrated over 
some (p t , 77) acceptance, a sum over Dirac delta functions 
(particle positions) 

n 

p(0)=5> l( 5(0-0,;) (1) 
z=l 



The Ti are weights (1, pt or E t ) appropriate to a given 
physical context. We assume integration over one unit of 
pseudorapidity, so multiplicity n estimates dn/dij (simi- 
larly for p t and E t ). The continuum parent density, not 
directly observable, is the object of analysis. Fixed parts 
of the parent density are estimated by a histogram aver- 
aged over an event ensemble. The discrete nature of the 
sample distribution and its statistical character present 
analysis challenges which are one theme of this paper. 

The correlation structure of the single-particle azimuth 
density is manifested in ensemble-averaged multiparticle 
(two-particle, etc.) densities. Accessing that structure 
by projection of multiparticle spaces to 2D or ID with 
minimal distortion is the object of correlation analysis. 
In each event the two-particle density is the Cartesian 
product p{4>i,4>2) = p(<pi) p{<h) 

n 

p(<h.,<h) = ( 2 ) 

2=1 

n.n—1 

+ ^ r iTj5{(j)l- (j>i)6{())2- (j>j), 
&j 

where the first term represents self pairs. Ensemble- 
averaged two-particle distribution p(4>i, 4>2) with correla- 
tions is not generally factorizable. By comparing the av- 
eraged two-particle distribution to a factorized (or mixed- 
pair) statistical reference two-particle correlations are re- 
vealed. In a later section we compare multipole moments 
from two-particle correlation analysis on azimuth to re- 
sults from conventional flow analysis methods. 



B. Fourier transforms on azimuth 

A Fourier series is an efficient representation if azimuth 
structure approximates a constant plus a few sinusoids 
whose wavelengths are integral fractions of 27r. A Fourier 
representation of a peaked distribution (e.g., jet cone) 
is not an efficient representation. We assume a simple 
combination of the lowest few Fourier terms. 

The Fourier forward transform (FT) is 

M) = z] — ^exp(im0) (3) 

m= — oo 

= + 2 coB(m[* - *„,]), 

m=l 

where boldface Q m is an event-wise complex amplitude, 
Q m is its magnitude and ^ m its phase angle. The second 
line arises because p(<p) is a real function, and the prac- 
tical upper limit on m is particle number n (wavelength 
~ mean interparticle spacing). Q m /2n = p m is the am- 
plitude of the density variation associated with the m th 
sinusoid. With — > 1 Q m is the corresponding num- 
ber of particles in 2ir if that density were uniform. ^ m 
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is event-wise by definition and does not require a tilde. 
The reverse transform (RT) is 



(4) 



Q m = / d(f>p(<f))exjp(-im<i 

J —IT 

n 

= y~]n exp(-zm<fc) 

i 

= Q m exp(im$ m ). 



For the discrete transform <f) £ [— 7r, 7r] is partitioned 
into M equal bins with bin contents fj and bin centers 
</);. We multiply Eq. © by bin width 5<f) = 2n/M, and 
the Fourier transform pair becomes 



n 



n M/2 n 

2 — cos ( w t 



,1/ 



71/ 



>j - * m ]) (5) 



m— 1 
M/2 

^ fj exp(-im0i), 

J=-M/2 



where the f; are also random variables. The upper limit 
M/2 in the first line is a manifestation of the Nyquist 
sampling theorem [1 71 ] . With r, — > f n — > and Q m /M 
is the maximum number of particles in a bin associated 
with the m th sinusoid. 



C. Autocorrelations and power spectra 



The azimuth autocorrelation density pa (<M) is a pro- 
jection by averaging of the pair density on two-particle 
azimuth space (fa, fa) onto difference axis fa\ — fa — fa. 
The autocorrelation concept is not restricted to peri- 
odic or bounded distributions or discrete Fourier trans- 
forms [lg]. In what follows pa includes self pairs. The 
autocorrelation density is defined as [l8j 



Pa(<PA) 



2tt 



dfao(4>) p(<j> + fa\) 



— X! Vir 3 \ 4>i)5(fa- 
1 ™ 

r rir i s ^ ~ & + ^ A )- 



(6) 



fa + fa\) 



Using Eq. (JSJ) we obtain the FT as 
1 f 7 " 

Pa(<Pa) = — J d(f> x 

X! ^ Zi exp(im0)x 



(7) 



m— — oo 



£ SkexpHm'[0 + ,M]) 



m' — — oo 



m— — oo 

r) 2 00 r> 2 



and the RT as 

Q m = n 2 (rcos(m[fa-^ m })) 2 



(8) 



= 2tt 



#a P~A(fa\) COs(7710a) 



= ^ nrj cos(m[0i - fa]) 

i,3'=l 

= rt(r 2 ) + n(n — l)(r 2 cos(m(/>A))- 

= n(r 2 ) + n(n - l)(r 2 cos 2 (m[0 - * r ])). 

Phase angle has been eliminated, and \I/ r will be 
identified with the reaction plane angle. 

We can write the same relations for ensemble-averaged 
quantities because the terms are positive-definite 



PaIPa) 
with RT 

Qm = 



Ql 



■ cos(m<pi 



2n 



[2^] 2 ^ [2tt] 

d(f>A Pa((^a) cos(m^A) 



(9) 



(10) 



nrj cos(m[0i - fa]) 



= n(r 2 ) + n(n — l)(r 2 cos(mfa\)) . 

We have adopted the convention Q 2 n = Q 2 ^ to lighten 
the notation. Coefficients Q m are power- spectrum ele- 
ments on wave-number index m. That FT transform pair 
expresses the Wiener- Khintchine theorem which relates 
power-spectrum elements Q m to autocorrelation density 
PA (0a)- The autocorrelation provides precise access to 
two-particle correlations given enough collision events, no 
matter how small the event-wise multiplicities. 



D. Autocorrelation structure 

The autocorrelation concept was developed in response 
to the Brownian motion problem and the Langevin equa- 
tion, a differential equation describing Brownian motion 
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which contains a stochastic term. The concept is al- 
ready apparent in Einstein's first paper on the subject 
(cf. App.jAj. The large-scale, possibly-deterministic mo- 
tion of the Brownian probe particle must be separated 
from its small-scale random motion due to thermal col- 
lisions with molecules. Similarly, we want to extract az- 
imuth correlation structure persisting in some sense over 
an event ensemble from event-wise random variations. 
The autocorrelation technique is designed for that pur- 
pose. 

The statistical reference for a power spectrum is the 
white-noise background representing an uncorrelated sys- 
tem. The reference is typically uniform up to large 
wave number or frequency (hence white noise), an in- 
evitable part of any power spectrum from a discrete pro- 
cess (point distribution). The "signal" is typically limited 
to a bounded region (signal bandwidth) of the spectrum 
at smaller frequencies or wave numbers. 

From Eq. ([5]) the power-spectrum elements are 

n n,n— 1 

Q 2 m = COs(m[& -<£,-]) (11) 

i i^j 

— n(r 2 ) + n(n — l)(r 2 cos(m0A)) 

The first term in Eq. (fTTj) is Q 2 e ^ the white- noise back- 
ground component of the power spectrum common to 
all spectrum elements. The second term, which we de- 
note V^, represents true two-particle azimuth correla- 
tions. Note that Q 2 , = n 2 (r 2 ), whereas V 2 = n(n—l)(r 2 }. 
In terms of complex (or vector) amplitudes we can write 

Q m = Qre/+V m , (12) 

where Q re / represents a random walker. There is no 
cross term in Eq. (|lip because Q re / and V m are uncor- 
related. 

Inserting the power-spectrum elements into Eq. §§§ we 
obtain the ensemble-averaged autocorrelation density 

PAW A) = -z d(0A) + ( 13 ) 

+ 2 E T53J cos (-^). 

m— 1 ^ ' 

The first term is the self-pair or statistical noise term, 
which can be excluded from pa by definition simply 
by excluding self pairs. The second term, with V 2 = 
n(n — l)(r 2 ), is a uniform component, and the third term 
is the sinusoidal correlation structure. The self-pair term 
is referred to in conventional flow analysis as the "auto- 
correlation," in the sense of a bias or systematic error, 
but that is a notional misuse of standard mathematical 
terminology. The true autocorrelation density is the en- 
tirety of Eq. (fT5)) , including (in this simplified case) the 
self-pair term, the uniform component and the sinusoidal 
two-particle correlations. 



In general, the single-particle ensemble-averaged dis- 
tribution po may be structured on (r],(j)). We want to 
subtract the corresponding reference structure from the 
two-particle distribution to isolate the true correlations. 
In what follows we assume 7"j = 1 for simplicity, therefore 
describing number correlations. We subtract factorized 
reference autocorrelation pA,ref{4>\,4>2) = Po(4>i) Po(<fe) 
representing a system with no correlations, with po = 
n/2n ~ d 2 n/dr)d4> in this simple example, to obtain the 
difference autocorrelation 

&Pa(<Pa) = PA~PA,ref (14) 
2 - 00 ta2 

= ^ + 2 E[^ c -(^a). 

m— 1 1 1 

The first term measures excess (non-Poisson) multiplicity 
fluctuations in the full (pt,r),(j)) acceptance. The second 
term is a sum over cylindrical multipoles. We now divide 
the autocorrelation difference by y/pA,ref = Po = n/2n 
to form the density ratio 

= 4^ +2V^ cos(m M (15) 

V m— 1 

A PA [0] , A PA [m] 
= + 2 X , cos ( m 

y/PA,ref ^ y/PA,ref 

The first term ApA [0]/ -JpA,ref is the density ratio av- 
eraged over acceptance (Arj,A<p). Its integral at full 
acceptance is normalized variance difference Aa 2 , 
(a^ — n) /n, which we divide by acceptance integral 1x2tt 
to obtain the mean of the 2D autocorrelation density 
Ap A [0]/^/pA~^f = Ao-l/iArj, A(f>)/AriA(p. The sinusoid 

amplitudes are ApA [m]/ \/pA/ref = n ( n — 1) v 2 n /(2Trn) = 
Po^mi defining unbiased v m . The event-wise = 
(cos(?ti0a)) are related to conventional flow measures 
v m , but may be numerically quite different for small 
multiplicities due to bias in the latter. Different mean- 
value definitions result in different measured quantities 
(cf. App. |B]). Whereas the Q 2 m are power-spectrum ele- 
ments which include the white-noise reference, the (oc 
squares of cylindrical multipole moments) represent the 
true azimuth correlation signal. This important result 
combines several measures of fluctuations and correla- 
tions within a comprehensive system. 

E. Azimuth vectors 

Power-spectrum elements Q 2 m are derived from com- 
plex Fourier amplitudes Q m . In an alternate representa- 
tion the complex Q m can be replaced by azimuth vectors 
Q m . The Q m are conventionally referred to as "flow vec- 
tors" [22J , but they include a statistical reference as well 
as a "flow" sinusoid. The V m defined in this paper are 
more properly termed "flow vectors," to the extent that 
such terminology is appropriate. We refer to the Q m by 
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the model-neutral term azimuth vector and define them 
by the following argument. 

The cosine of an angle difference — cos(m[(j>i — fa]) = 
cos(m^i) cos(m02) + sin(m0i) sin(m02) — can be repre- 
sented in two ways, with complex unit vectors u(m<j)) = 
exp(im(j)) or with real unit vectors u(m(f>) = cos(m</>)z + 
sm(m(f>)j [complex plane (5Rz,9z) vs real plane (x,y)]. 
Thus, 

cos(m[^i - fa]) = !ft{u(m</>i) u*(mfa)} (16) 
= u(mfa) ■ u(mfa)). 

If an analysis is reducible to terms in cos(rn[0i — fa]) 
the same results are obtained with either representation. 
Thus, we can rewrite the first line of Eq. ([3|) as 

ft® = E ^- 5 H) (it) 

m— — 00 

- g +2 f:^cosM0-* m ]), 

TCI — 1 

in which case 

d<p p((f>) u{m(j>) = nujmfa) (18) 

= n(r[cos(m0), sin(m</))]) 
= Q m u(m^ m ), 

an event- wise random real vector. 



III. CONVENTIONAL METHODS 

We now use the formalism in Sec. |TT] to review con- 
ventional flow analysis methods in a common framework. 
We consider five significant papers in chronological order. 
The measurement of angular correlations to detect col- 
lective dynamics (e.g., parton fragmentation and/or hy- 
drodynamic flow) proceeds from directivity (1983) at the 
Bevalac (1-2 GeV/u fixed target) to transverse spheric- 
ity predictions (1992) for the SPS/RHIC (^/JWn = 17 - 
200 GeV), then to Fourier analysis of azimuth distribu- 
tions (1994, 1998) and v 2 centrality trends (1999). We 
pay special attention to the manipulation of random vari- 
ables (RVs). RVs do not follow the algebra of ordinary 
variables, and the differences are especially important for 
small sample numbers (multiplicities, cf. App. IB"]) . 

A. Directivity at the Bevalac 

An important goal of the Bevalac HI program was ob- 
servation of the collective response of projectile nucleons 
to compression during nucleus-nucleus collisions, called 
directed flow, which might indicate system memory of 
the initial impact parameter as opposed to isotropic ther- 
mal emission. Because of finite (small) multiplicities and 



large fluctuations relative to the measured quantity the 
geometry of a given collision may be poorly defined, but 
the final state may still contain nontrivial collective infor- 
mation. The analysis goal becomes separating possible 
collective signals from statistical noise. 

An initial search for collective effects was based on the 
3D sphericity tensor S = Y^i PiPi 0, EH described in 
App. lC 31 Alternatively, the directivity vector was defined 
in the transverse plane [lj|. In the notation of Sec. [Til 
including weights — > Wipu, directivity is 

n n 
Qi = y~] v-kpti = ^2 Wip ti u(fa) (19) 

i i 

= QimO&i), 

azimuth vector Q m with m = 1 (corresponding to di- 
rected flow) . Event-plane (EP) angle estimates true 
reaction-plane (RP) angle To maintain correspon- 
dence with SPS and RHIC analysis we simplify the de- 
scription in [l9| to the case that n single nucleons are 
detected and no multi-nucleon clusters; thus a — > 1 and 
A-*n. 

The terms in Qi are weighted by u>i — w(yi) — ±1 cor- 
responding to the forward or backward hemisphere rela- 
tive to the CM rapidity, with a region [—6, S] about mid- 
rapidity excluded from the sum (wi = 0). Qi then ap- 
proximates quadrupole moment q2i derived from spher- 
ical harmonic SiKj ^ sin(2#) cos(0), as illustrated in 
Fig. [TJ (left panel, dashed lines). In effect, a rotated 
quadrupole is modeled by two opposed dipoles, point- 
symmetric about the CM in the reaction plane. It is ini- 
tially assumed that EP angle ^1 of vector Qi estimates 
RP angle , 3> r , and magnitude Qi measures directed flow. 

The first part of the analysis was based on subevents — 
nominally equivalent but independent parts of each 
event. The dot product Qia ■ Qib for subevents A 
and B of each event was used to establish the exis- 
tence of a significant flow phenomenon and the angu- 
lar resolution of the RP estimation via the distribu- 
tion on — ^ib- The EP resolution was defined by 

cos(>3>i - * r ) = 2y / cos(* 1A - ^>ib)- The magnitude of 

Qi was then related to an estimate of the mean trans- 
verse momentum in the RP. Integrating over rapidity 
with weights Wi we obtain event-wise quantities 

n n,n—l 
Ql = E W i P ^ + E Wl W 3 P U ' ( 20 ) 

~ n(p 2 t ) + n{n - l)(p 2 t cos(<£a)) 

= Qlef + V?. 

The last line makes the correspondence with the notation 
of this paper. 

The initial analysis in [lj| used Q\ — Q 2 ef = V 2 = 
n(n— l)(p 2 ), assuming that x is contained in the RP and 
there are no non-flow correlations. Note that n(n — 1) = 



7 



ESy 1 | u, » u 'jl contains weights WiWj implicitly. Since 
Wi ~ sin[2 0(yj)] (cf. Fig. Q]- left panel), what is actu- 
ally calculated in [19( for the single-particle (no multi- 
nucleon clusters) case is the pt-weighted r.m.s. mean 
of the spherical harmonic SRY^ 0) quadrupole mo- 
ment q2i, thereby connecting rank-1 tensor Qi (with 
weights on rapidity) to rank- two sphericity tensor 5 (cf. 
App. [C"3j) . The mean-square quantity calculated is 



Pi 



sin 2 (26»)p t 2 cos 2 (0- * r ) 



n(n — 1) 



sin 2 (26>) 



(21) 



a minimally-biased statistical measure as discussed in 
App. [Bj from which we obtain p x — \fp\ estimating the 
transverse momentum per particle in the reaction plane. 

The second part of the analysis sought to obtain p x (y) , 
the weighted-mean transverse momentum in the RP as a 
function of rapidity. It was decided to determine p X i — 
Pu cos(0i — ^v) for the i th particle relative to the RP, but 
with \fv estimated by EP angle The initial attempt 
was based on 



w l p' x . 



WiPti ■ = WiPti 



En -> 
j W 3 Ptj 



T,k w kPtk\ 

Summing Wi pu over all particles in a y bin gives 

V^n 2 2 y,n,n-' 
Z^i w i Pti + l^ijtj 



(p' x ) = 



WiWjPti -p t j 



E" M IE* wfcPtfel 

n(p 2 ) + V? 
nQi 

Qi/n = 



(22) 



(23) 



(P?) 



1 



from which we obtain ensemble mean p^, = (pf,.). That 
result can be compared directly with the linear speed in- 
ferred from a random walker trajectory, which is 'infinite' 
in the limit of zero time interval (cf. App. |A"|) . The first 
term in the radicand is said to produce "multiplicity dis- 
tortions" in conventional flow terminology. The second 
term contains the unbiased quantity. 

In contrast to the first part of the analysis the sec- 
ond method retains the statistical reference within Qi 
as part of the result, so that Qi/n ~ y (pf) /n for small 
multiplicities and/or flow magnitudes, a false signal com- 
parable in magnitude to the true flow signal. The unbi- 
ased directed flow p x was said to be "distorted" by the 
presence of the statistical reference (called unwanted self- 
correlations or "autocorrelations" ) to the strongly-biased 
value p' x dominated by the statistical reference. 

An attempt was made to remove statistical distortions 

arising from self pairs by redefining Qi — ► Qu, a vector 
complementary to each particle i with that particle omit- 
ted from the sum. The estimator of \I/ r for particle i is 

then *h in Qu = YJ^l w jPtj = Qn and 



WiPxi 



En — 1 -* 
&i w iPtj 

\TZ^l W kPtk\ 



.(24) 
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FIG. 1: Left panel: Comparison of directed flow data from 
Fig. 3(a) of the event-plane analysis and the 5ft Y2 1 °c 
sin(29[yi a b]) spherical harmonic, with amplitude 95 MeV/c 
obtained from V\ . Weights in the form w(yi a b) are denoted 
by dashed lines. The correspondence with sin(28[yi a t]) (solid 
curve) is apparent. Right panel: The EP resolution obtained 
from [22I ] (dashed curve) and from ratio \J (n —- l)V 1 2 /nQ' 1 2 
defined in this paper (solid and dotted curves for several n 
values) . 



Summing over i within a rapidity bin one has 

n \-^n — 1 -* -* 

1 x - Ej^i " <"'.;/'•<< •/'••.. 



En^ 

XL 



Efc^i w k ptk\ 



(25) 



= V>!> 



^/(pi) + (n-2)(pl) 



l)(p 2 ) + (n-l)(n-2)(p 2 cos(0 A )>. 



with Qi = \J (n 
Since — n(n—l)(p x ) = n(n — l)(p? cos(</>a)) one 
sees that the division by Q\ is incorrect, even though 
it seems to follow the chain of argument based on RP 
estimation and EP resolution with correction. The 

correct (minimally-biased) quantity is p x = 



, pI = 

y Vi /n(n — 1). The new EP definition removes the ref- 
erence term from the numerator, but Q[ in the denom- 
inator retains the statistical reference in p". There is 

the additional issue that Vx 2 ^ x. Two different mean 

^ The dif- 



values are represented by (p x ) and p x = \i p x 
ference can be large for small event multiplicities. 

The remaining bias was attributed to the EP resolu- 
tion. The resolution correction factor derived from the 
initial subevent analysis (cf. App. [Dj was applied to 
(p") to further reduce bias. In Fig. [1] (right panel) we 
compare the EP resolution correction from [22| (dashed 



curve) with factor y (n — l)Vj /nQf required to convert 

(p' x ) from Eq. ([25]) to p x from Eq. (|2ip. The agreement 
is very good. Eq. ([2Tjl is the least biased and most direct 

way to obtain p x 



ip x , both globally over the detec- 
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tor acceptance and locally in rapidity bins, without EP 
determination or resolution corrections. 

In Fig. Q] (left panel) we show the data (points) for 
Px(yiab) from the EP-corrected analysis and the solid 
c urve p x sm[29 (yi ab )} oc q 2 iY 21 (9(y), 0), where p x = 

\Jv 2 /n(n — 1) = 95 MeV/c. The agreement is good, 
and the similarity of sin[20(y; a f,)] (solid curve) to weights 
w(yiab) (dashed lines) noted above is apparent. Loca- 
tion of the sin[20(y/ a (,)] extrema near the kinematic limits 
(vertical lines) is an accident of the collision energy and 
nucleon mass. These results are for E cm = 1.32 GeV ~ 
In the notation of this paper V 2 — 4.7 

0.095 GeV/c 



\f~2mNC 2 . 

(GeV/c) 2 , yjv?/n(n-l) 



wpx/a, and 



Q x = n\JV?ln(n-l) = 2.17 GeV/c ~ ( not 
By direct and indirect means (directivity and RP esti- 
mation) quadrupole moment q2i oc V\ was measured. 

B. Transverse sphericity at higher energies 

The arguments and techniques in (2(| suggest a smooth 
transition from relativistic Bevalac and AGS energies 
(collective nucleon and resonance flow) to intermediate 
SPS and ultra-relativistic RHIC energies (possible trans- 
verse flow, possibly QCD matter, collectivity manifested 
by correlations of produced hadrons, mainly pions). For 
all heavy ion collisions thermalization is a key issue. 
Clear evidence of thermalization is sought, and collective 
flow is expected to provide that evidence. 

Two limiting cases are presented in [20] for SPS flow 
measurements: 1) linear N-N superposition with no col- 
lective behavior (no flow); 2) thermal equilibrium - col- 
lective pressure in the reaction plane - fluid dynamics 
leading to "elliptic" flow. In a hydro scenario the initial 
space eccentricity transitions to momentum eccentricity 
through thermalization and early pressure. The paper 
considers flow measurement techniques appropriate for 
the SPS and RHIC, and in particular newly defines trans- 
verse sphericity St- 

According to 12011 the 3D sphericity tensor introduced 
at lower energies |15| can be simplified in ultra-relativistic 
heavy ion collisions to a 2D transverse sphericity tensor. 
Sphericity is transformed to 2D by p — ► pt , omitting the 
momentum z component near mid-rapidity. Transverse 
sphericity (in dyadic notation) is 



2S t = 2 pupti 

i 

i 

= n(p 2 t ){l + a 1 C{^ 2 )} 



(26) 



defining a% and ^ 2 in the tensor context, with 

c(4>) 



cos(2</>) sin(2^>) 
sin(20) — cos(2c/>) 



This a definition corresponds to Eq. (3.1) of 
We next form the contraction of St with itself 

2S t :S t = 2j2(Pu-p tl ) 2 



(28) 



= 2 Pu + 2 cos2 - h ) 

i i£j 
= 2n(pf) + 2n(n - l)(p* cos 2 (0 A )) 
= n(n + l)(p*) + n(n - l){p 4 t cos(20 A )) 

using the dyadic contraction notation A : B = AabBab, 
with the usual summation convention. That self- 
contraction of a rank-2 tensor can be compared to the 

more familiar self-contraction of a rank-1 tensor Q 2 -Q 2 = 
Q\ = n(p 2 ) + n(n — l)(p 2 cos(20a))- The quantity 
2[St : St] re f = n(n + l)(Pt) is the (uncorrelated) refer- 
ence for the rank-2 contraction, whereas Q\ re j — n(p 2 ) 
is the reference for the rank-1 contraction. Subtracting 
the rank-2 reference contraction gives 

2 {<St : S t - [S t : St] re f \ = n(n - l)(p£ cos(20a)) 

* (Pt)V 2 2 , (29) 

which relates transverse sphericity to two-particle corre- 
lations in the form V 2 = Q\ — Q\ re p thus establishing 

the exact correspondence between St and Q 2 . 

From the definition of a\ in Eq.[5S]above and Eq. (3.1) 
of [201 we also have 



2Sf.S t = n 2 {p 2 ) 2 {l + a\) 



which implies 



n 2 {p 2 ) 2 &\ = n 2 a 2 p2t +n(pt) 

+ n(n - l)(pt cos(2^ A )) 
- n 2 d 2 p2t +{p 2 )Q 2 2 . 



(30) 



(31) 



(32) 



The first relation is exact, given the definition of a%, 
but produces a complex statistical object containing the 
event-wise variance of p\ in its numerator and random 
variable n 2 in its denominator. For X/n — > it is true 
that 5i\ — > (cos(2[0— ^ r ])) (since ^ 2 — > ^ r also), but 6>x 
is a strongly biased statistic for finite n. 
The definition 



«2 



(Pi ~ P\) 

(pI+pI) 



(33) 



(27) 



from Eq. (2.5) of [20| seems to imply a 2 — (p 2 cos(2[0 — 
^r]))/(p 2 } — > (cos(2[</>— ^t-])), assuming that x lies in the 
RP. However, the latter relation fails for finite multiplic- 
ity [the effect of the statistical reference or self-pair term 
n in Eq. f3T]) ] because each of event-wise (p 2 ) and (p 2 ) is 
a random variable, and their independent random varia- 
tions do not cancel in the numerator. The exact relation 
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is the first line of 

n 2 (p 2 ) 2 a 2 = n 2 (p 2 t cos^-^})) 2 (34) 
~ n 2 (p 2 )(p t cos(2[cj ) ^^ 2 })) 2 

The second line is an approximation which indicates that 
«2 is more directly related to Q2 than is c\\. But Q\ is a 
poor substitute for V 2 2 which represents true two-particle 
azimuth correlations in a minimally-biased way by incor- 
porating a proper statistical reference. The effect of the 
reference contribution is termed a 'distortion' in 201. 



C. Fourier series I 

Application of Fourier series to azimuth particle dis- 
tributions was introduced in (2lj . Fourier analysis is de- 
scribed as model independent, providing variables which 
are "easy to work with and have clear physical interpreta- 
tions." Sinusoids or harmonics are associated with trans- 
verse collective flow, the model-dependent language fol- 
lowing [2(| closely. To facilitate comparisons we convert 
notation in to that used in this paper: r((f>) — > p(4>), 

(^m?ym) * Qm: ^Pm * ^ mi ^m * Qm and V m ► V m . 

According to the proposed method density p(4>) repre- 
sents, within some (pt, n) acceptance, an event- wise par- 
ticle distribution on azimuth (f> including weights n = 
1, pu or E t i- The FT in terms of azimuth vectors is 

n 

P{4>) = J2 TiS ^-^ ( 35 ) 

i 

rn—l 
m=l 

and the RT is 

_ n 

Qm = S^nujmfa) = Q m u(m^ m ), (36) 

i 

forming a conventional Fourier transform pair [cf. 
Eqs. ^ and ||3J)]. Scalar amplitude Q m = r iu(m(j>i) ■ 
ulm^m) = n(r cos(to[</> — 'I'm])) is the proposed flow- 
analysis quantity. Q m is said to measure the flow mag- 
nitude, and \& m estimates reaction-plane angle f r . It 
is proposed that Q m (ii) evaluated within bins on n may 
characterize complex "event shapes" [densities on (77, <fi)]. 

As with directivity and sphericity, multiplicity fluctu- 
ations are seen as a major obstacle to flow analysis with 
Fourier series. Finite multiplicity is described as a source 
of 'bias' which must be suppressed. It is stated that 
(V m , *& r ) are the "parameter^] relevant to the magnitude 



of flow," whereas the observed (Q m , ^m) are biased flow 
estimators. In the limit 1/n — » the two cases would 
be identical. A requirement is therefore placed on min- 
imum event-wise multiplicity in a "rapidity slice." To 
solve the finite-number problem the paper proposes to 
use the event frequency distribution on Q 2 ^ from event- 
wise Fourier analysis to measure flow. 
If correlations are zero then 

n 

Q 2 m ^Ql e! = ^=n{r 2 )^a\ (37) 

i 

with (r 2 ) = <Tq. By fitting the frequency distribution on 
with a model function it is proposed to obtain V m as 
the unbiased flow estimator. The fitting procedure is said 
to require sufficiently large event multiplicities to obtain 
V m unambiguously. 

The distribution on Q 2 n is derived as follows (cf. 
Fig. [2] - left panel). The magnitude of statistical ref- 
erence Qref (a random walker) has probability distri- 
bution oc exp(-Q 2 e/ /2(3 2 e/ ), with Q 2 e/ = n(r 2 ). But 

Qref = Qm — V m , therefore (cf. Fig. [5]- left panel) 

Qlef = Qrn + V m - 2Q m ^ m cos(m[tf m - r ]), (38) 

and 

exp(-Q 2 e/ /2Q 2 e/ ) p(Q m , tf m ; V m ,V r ). (39) 

When integrated over cos(m[^I' TO — \Dv]) there results the 
required probability distribution on Q 2 m , with fit param- 
eter V m - The distribution on Q 2 n is said to show a 'non- 
statistical' shape change from which V m can be inferred 
by a model fit "/ree from uncertainties in event-wise de- 
termination of the reaction plane." It is also proposed 
to use p(Qm,^m',V m ,^ r ) to determine the EP resolu- 
tion cos^m — \fv) by integrating over Q 2 m and using the 
resulting projection on cos(^ m — \fv) to determine the 
ensemble mean (Fig. 3 of [2l|). 

While one could extract ensemble-average V m at some 
level of accuracy by fitting the frequency distribution on 
Qm with a model function, we ask why go to that trou- 
ble when V m is easily obtained as a variance difference? 
Instead of Eq. (j38|) we simply write 

QL = Qlef + V 2 , (40) 

where the cross term is zero on average and Q 2 e j — n(r 2 ) 
represents the power-spectrum white noise, which is the 
same for all m (i.e., 'white'). If the vector mean values 
are zero that is a relation among variances. Ensemble 
mean V 2 n = Q 2 n — Q 2 re f is therefore simply determined. 

For the EP resolution we factor V m 

V 2 = n(n - l)(r 2 cos 2 (m[0 - * r ])) (41) 
= n(n - l)(r 2 cos 2 (rn[0 - tf m ])) cos 2 (m[# m - *,.]). 
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We use Q m — n(r cos(m[0 — W m ])) and the assumption 
that random variable ty m — \& r is uncorrelated with <fi — 
to obtain 



cos 2 (m[# m - * r ]) = 



(n-l)i 



(42) 



which defines the EP resolution of the full n-particle 
event in terms of power-spectrum elements (cf. App. |Pj) . 
The square root of that expression is plotted in Fig. 
(right panel) as the solid curves for several values of n. 
The solid and dotted curves are nearly identical to those 
in Fig. HJ 






1 




[ i | i i i i | i , w... 




0.9 






1 


0.8 






= 


0.7 








0.6 






B 


0.5 




/ "V„ i 


COS 


0.4 




/ V(n-l)Q ra j 


V 


0.3 




/ n-5 




0.2 




n = 10,... ,50 j 




0.1 




[22] 1 














v /a ■■ 



V(2nv2)-V(2V» 3 



FIG. 2: Left panel: Distribution of event-wise elements of Q m 
components determined by the gaussian-distributed random 

walker Q re f and possible correlation component V m . Right 
panel: Reaction-plane resolution estimator (cos(m<5\E , mr )), 
with S^mr = — ^ r , determined from fits to a distribution 
on Q„ as in the left panel (dashed curve) , and from Eq. [42] 
for several values of n (solid and dotted curves). 



As in other flow papers there is much emphasis on 
insuring adequate multiplicities to reduce bias to a man- 
agable level, because an easily-determined statistical ref- 
erence is not properly subtracted to reveal the contribu- 
tion from true two-particle correlations in isolation. For 
V2nv m > 1 in Fig. [2] (right panel) the EP is meaningful; 
bias of event-wise quantities relative to the EP is man- 
ageable. For \f2nv m < 1/2 the EP is poorly defined, 
and ensemble-averaged two-particle correlations are the 
only reliable measure of azimuth correlations. In cither 
case EP estimation is only justified when a non-flow phe- 
nomenon is to be studied relative to the reaction plane. 



D. Fourier series II 

A more elaborate review of flow analysis methods 
based on Fourier series is presented in [22J. The ap- 
proach is said to be general. The event plane is ob- 
tained for each event. The event-wise Fourier ampli- 
tude^) q m = Q m /Qo relative to the EP are corrected for 
the EP resolution as obtained from subevents. We mod- 
ify the notation of the paper to v!j£ s — ► q m and to, — > fj to 
maintain consistency within this paper. We distinguish 
between the unbiased v m = V m /Vo and the biased q m . 



According to |22|], in the 1/n — ► limit (Q m — > V m , no 
tildes, no random variables) the dependence of the single- 
particle density on azimuth angle <j> integrated over some 
(p t , y) acceptance can be expressed as a Fourier series of 
the form 



2 J2 v r, 



cos m 



0-* r )] 



(43) 



with reaction-plane angle ^f r . As we have seen, the factor 
2 comes from the symmetry on index m for a real-number 
density, not an arbitrary choice as suggested in (22J. In 
this definition Vq/2it has been factored from the Fourier 
series in [2l|. The Fourier "coefficients" v m in this form 
(actually coefficient ratios) are not easily related to the 
power spectrum. In the 1/n — ► limit the coefficients 
are v m = (cos(m[(j) — \fv]))- 

In the analysis of finite-multiplicity events, reaction- 
plane angle defined by the collision (beam) axis and 
the collision impact parameter is estimated by event 
plane (EP) angle $> m , with $ m derived from event- wise 

azimuth vector Q m (conventional flow vector) 

~ n 

Qm = y^jSimfa) = Q m u(m^f m ). (44) 
The finite-multiplicity event-wise FT is 



2 9» 



COS 777 



(45) 



with q m = (cos(m[</>-* m ])>; e.g., q 2 ~ a 2 (cf. Eq. (JM])). 

According to the conventional description the EP an- 
gle is biased by the presence of self pairs, unfortunately 
termed the "autocorrelation effect" or simply "autocor- 
relation" in conventional flow analysis p^ . [23 | , whereas 
autocorrelations and cross-correlations are distributions 
on difference variables used for decades in statistical anal- 
ysis to measure correlation structure on time and space. 
As in Jj|, to eliminate "autocorrelations" EP angle 
is estimated for each particle i from complementary flow 
vector Q mi = YTj^l r i^( 2 ^j) = Qmiu(m^ m i), a form of 
subevent analysis with one particle vs n — 1 particles (cf. 
App.0. 

In the conventional description the event-plane res- 
olution results from fluctuations 5^> r = 9 m — ^> r of 
the event-plane angle fy m (or $ m i) relative to the true 
reaction-plane angle VEV (e.g., due to finite particle num- 
ber). The EP resolution is said to reduce the observed 
q m relative to the true value v m : 



1 



1 



77 77 

cos (m[^ r 



(46) 



The EP resolution (first factor, second line) is obtained in 
a conventional flow analysis in two ways: the frequency 
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distribution on $ m — fy r discussed in Sec. IIII CI and the 
subevent method discussed in App. [D] A parameteri- 
zation from the frequency-distribution method reported 
in [22j is plotted as the dashed curve in Fig. [5] (right 
panel) . There is good agreement with the simple expres- 
sion y/n/ (n — 1) V m /Q m obtained in Eq. (|4"2"|) , which also 
follows from Eq. (|4*6")l . 

Two methods are described for obtaining v m without 
an event-plane estimate, with the proviso that large event 
multiplicities in the acceptance are required. The first, in 
terms of the conventional flow vector, is expressed (with 
n -> 1) as (Eq. (26) of (H) 

Ql = n + rfv 2 m . (47) 

That expression is constrasted with the exact event-wise 
treatment, where for each event we can write 

Q 2 m = n + 53cos[m(&-&)] ( 48 ) 

= n + n(n — 1)(cos(7ti0a)) 
= n + n(n- l)v^ 
-> Qm = n + n(n - 1)5^ = n + V 2 

Note the similarity with fluctuations measured by num- 
ber variance <r^ = u+A(T„, where the second term on the 
RHS is an integral over two-particle number correlations 
and the first is the uncorrelated Poisson reference (again, 
the self-pair term in the autocorrelation). 

There are substantial differences between the two 
formulations above, especially for smaller multiplicities. 
Eq. (|48[l is unbaised for all n and provides a simple way 
to obtain = n{n — 1)5^, = Q 2 m — n. The conventional 
method uses a complex fit to the frequency distribution 
on Q 2 n to estimate V m as in [2l[ . Why do that when such 
a simple alternative is available? 

E. V2 centrality dependence 

In [23| the expected trend of v 2 with A-A centrality for 
different collision systems is discussed in a hydro context. 
It is stated that the v 2 centrality trend should reveal the 
degree of equilibration in A-A collsions. The centrality 
dependence of v 2 /e should be sensitive to the "physics of 
the collision" — the nature of the constituents (hadrons 
or partons) and their degree of thermalization or collec- 
tivity. "It is understood that such a state requires (at 
least local) thermalization of the system brought about 
by many rescatterings per particle during the system 
evolution... V2 is an indicator of the degree of equilibra- 
tion." 

Thermalization is related to the number of rescatter- 
ings, which also strongly affects elliptic flow according to 
this hydro interpretation. In the full hydro limit, corre- 
sponding to full thermalization where the mean free path 
A is much smaller than the flowing system, relation v 2 oc e 
is predicted, with e the space eccentricity of the initial 



A-A overlap region 20 . Conversely, in the low-density 
limit (LDL), where A is comparable to or larger than 
the system size, a different model predicts the relation 
v 2 oc eA 1 ' s /A, where A 1 / 3 /A estimates the mean num- 
ber of collisions per "particle" during system evolution 
to kinetic decoupling [2J|. In the LDL case v 2 oc £-5^, 
where S = irR x R y is the (weighted) cross-section area of 
the collision and e = R \ +R l is the spatial eccentricity. 
Those trends are further discussed in App. [E] 

According to the combined scenario, comparison of 
the centrality dependence of v 2 at energies from AGS 
to RHIC may reveal a transition from hadronic to ther- 
malized partonic matter. The key expectation is that at 
some combination(s) of energy and centrality v 2 /e tran- 
sitions from an LDL trend (monotonic increase) to hydro 
(saturation), indicating (partonic) equilibration. 

However, it is important to note two things: 1) That 
overall description is contingent on the strict hydro sce- 
nario. If the quadrupole component of azimuth correla- 
tions arises from some other mechanism then the descrip- 
tions in [2(| and are invalid, and v 2 does not reveal 
the degree of thermalization. 2) v 2 is a model-dependent 
and statistically-biased quantity motivated by the hydro 
scenario itself. The model-independent measure of az- 
imuth quadrupole structure is V 2 /n = nuf (defining an 
unbiased v 2 )- It is important then to reconsider the az- 
imuth quadrupole centrality and energy trends revealed 
by that measure to determine whether a hydro interpre- 
tation is a) required by or even b) permitted by data. 



IV. IS THE EVENT PLANE NECESSARY? 

A key element of conventional flow analysis is estima- 
tion of the reaction plane and the resolution of the esti- 
mate. We stated above that determination of the event 
plane is irrelevant if averaged quantities are extracted 
from an ensemble of event-wise estimates. The reaction- 
plane angle is relevant only for study of nonflow (minijet) 
structure relative to the reaction plane on fe, the sum 
(pair mean azimuth) axis of ((j>i,(j) 2 ). In contrast to con- 
cerns about low multiplicities in conventional flow analy- 
sis, proper autocorrelation techniques accurately reveal 
"flow" correlations (sinusoids) and any other azimuth 
correlations, even in p-p collisions and even within small 
kinematic bins. In this section we examine the necessity 
of the event plane in more detail. 

The reaction plane (RP), nominally defined by the 
beam axis and the impact parameter between centers of 
colliding nuclei, is determined statistically in each event 
by the distribution of participants. The RP is estimated 
by the event plane (EP), defined statistically in each 
event by the azimuth distribution of final-state particles 
in some acceptance. We now consider how to extract 
v m relative to a reaction plane estimated by an event 
plane in each event. Several different flow measures are 
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implicitely defined in conventional flow analysis 



v m = cos(m[0 — f r ]) ideal case, (49) 

= (cos 2 (m[(j) — \? r ])) unbiased estimate 

q m = (cos(m[(j) ~ ^ m ])) self-pair bias 

v' m = (cos(m[(j)i — $ m J)) reduced bias 

v' m is the event- wise result of a conventional flow analysis. 
The ensemble average v' m must be corrected for the "EP 
resolution" which we now determine. 

The basic event-wise quantities, starting with an inte- 
gral over two-particle azimuth space, are 



V, 



n,n—l 

" ~ u(m4>i) ■ u(m4>j) 



(50) 



= n(n — 1)(cos(t7i</>a)) 
= n(n — l)(cos 2 (m 



*r])> 



n(n - l)vl 



For the limiting case of subevents A and B with A a single 
particle the subevent azimuth vector complementary to 
particle i is 



The right panel indicates that for large nv 2 ^ single- 
particle reaction-plane estimates can provide a "flow" 
measurement with manageable bias. For small ni^ the 
EP resolution averaged over many events is itself a "flow" 
measurement, even though the reaction plane is inacces- 
sible in any one event. is determined from the same 
underlying two-particle correlations by other means — the 
only difference is how pairs are grouped across events. In 
App. |D] the EP resolution is determined for the case of 
equal subevents A and B 

From this exercise we conclude that event-plane deter- 
mination is irrelevant for the measurement of cylindrical 
multipole moments ("flows"). Following a sequence of 
analysis steps in the conventional approach based on de- 
termination of and correction for the EP estimate, the 
event plane cancels out of the flow measurement. What 
results from the conventional method is approximations 
to signal components of power-spectrum elements which 
can be determined directly in the form V m jn = nw„ 
obtained with the autocorrelation method, which defines 
an unbiased version of v m . Event-plane determination 
can be useful for study of other event-wise phenomena in 
relation to azimuth multipoles. 



Qmi = / ] = QmSim^mi). (51) 

We make the following rearrangment 

n n—l 

Vl = ]T i?(m&) • u(m^) (52) 

i j^i 



^ Qmi 



co$(m[<f>i - * mi ]) 



nQ' m (cos(m[((b-^ m })) 



n(n - l)v 2 m = nQ' m v' n 



n-1 V m 



and, since Q' m ~ sjn — 1 + (n — l)(n — 2)(cos(m0A)), 
we identify the EP resolution as 



(cos(m[*^ - * r ])) = 



(53) 



Qm 

V' 

m 
Q'm 



where the primes refer to a subevent with n—l particles. 
That expression, for full events with multiplicity n, is 
plotted in Fig. [2] (right panel) for several choices of n. An 
n-independent universal curve on V m fn is multiplied by 
n-dependent factor \/ n/(n — 1), where n is the number 
of samples in the event or subevent. The dashed curve is 
the parameterization from [221 ] . 



V. 2D (JOINT) AUTOCORRELATIONS 

We now return to the more general problem of angular 
correlations on (^1,772,^1,^2). We consider the analysis 
of azimuth correlations based on autocorrelations, power 
spectra and cylindrical multipoles without respect to an 
event plane in the context of the general Fourier trans- 
form algebra presented in Sec. ITT] We seek a comprehen- 
sive method which treats n and (f> equivalently. 

In conventional flow analysis there are two concerns 
beyond the measurement of flow in a fixed angular accep- 
tance: a) study flow phenomena in multiple narrow ra- 
pidity bins to characterize the overall "three-dimensional 
event shape," analogous to the sphericity ellipsoid but 
admitting of more complex shapes over some rapidity in- 
terval, and b) remove nonflow contributions to flow mea- 
surements as a systematic error. Maintaining adequate 
bin multiplicities to avoid bias is strongly emphasized in 
connection with a). The contrast between such individ- 
ual Fourier decompositions on single-particle azimuth in 
single-particle rapidity bins and a comprehensive analy- 
sis in terms of two-particle joint angular autocorrelations 
is the subject of this section 



A. Stationarity condition 

In Fig. [3] we show two-particle pair-density ratios f = 
p/pref on (771,772) (left panel) and (cj>i, (fe) (right panel) 
for mid-central 130 GeV Au-Au collisions [l^- The hat 
on f indicates that the number of mixed pairs in p re f 
has been normalized to the number of sibling pairs in 
p. In each case we observe approximate invariance along 



13 



sum axes rjs, = f]i + 772 and 4>t, — 4>i + 4>2- In time-series 
analysis the equivalent invariance of correlation structure 
on the mean time is referred to as stationarity, implying 
that averaging pair densities along the sum axes loses no 
information. The resulting averages are autocorrelation 
distributions on difference axes 77 a and 0a • 




FIG. 3: Normalized like-sign pair-number ratios f — p/p re f 
from central Au-Au collisions at 130 GeV for (771 , T72) (left 
panel) and (4>i,4>2) (right panel) showing stationarity — 
approximate invariance along sum diagonal xi + x%. 




FIG. 4: Autocorrelation averaging schemes on xt, = xi + X2 
for a prebinned single-particle space (left panel) and for pairs 
accumulated directly into bins on the two-particle difference 
axis (right panel). 



a given 2 a is the bin sum along that diagonal divided by 
the averaging interval. 



In Fig. [3] (right panel) one can clearly see the cos(20a) 
structure conventionally associated with elliptic flow 
(quadrupole component). However, there are other con- 
tributions to the angular correlations which should be 
distinguished from multipole components, accomplished 
accurately by combining the two angular correlations into 
one joint angular autocorrelation. 

B. Joint autocorrelation definition 

With the autocorrelation technique the dimensional- 
ity of pair density p(rji, 772, 0i, 02) can be reduced from 
4D to 2D without information loss provided the dis- 
tribution exhibits stationarity. Expressing pair density 
p{i)\ifh,<hi<h) K'7B) } 7Ai^s,0a) m differential form 
d 4 n/dx 4 we define the joint autocorrelation on (77A, 0a) 

I A \ - d2 / d2 "( ? 7s,?7A,0s,0A)\ 
Pa(va,M = — -7— { ^—7- > (54) 

by averaging the 4D density over t]y, and 0£ within a 
detector acceptance. The autocorrelation averaging on 
77s is equivalent to (dn/drj) v w n(Arj)/Ar/ at 77 = 0. The 
magnitude is still the 4D density d A n / '<i?7i<i?72(i0i<i02, but 
it varies only on the two difference axes (77 a , 0A ) ■ 

In Fig. H] we illustrate two averaging schemes [TH . In 
the left panel we show the averaging procedure applied 
to histograms on (x\,X2) as in Fig. [3] Index k denotes 
the position of an averaging diagonal on the difference 
axis. In the right panel we show a definition involving 
pair cuts applied on the difference axes. Pairs are his- 
togrammed directly onto ( 77a, 0a)- Periodicity on im- 
plies that the averaging interval on 0s /2 is 27r indepen- 
dent of 0a- However, 77 is not periodic and the averaging 
interval on rjs:/2 is A77 — |t7a|j where A77 is the single- 
particle rj acceptance [181 ] . The autocorrelation value for 



C. Autocorrelation examples 

In Fig. [S] we show p t joint angular autocorrelations for 
Hijing central Au-Au collisions at 200 GeV [Hj]. The left 
panel shows quench-on collisions. The right panel shows 
quench-off collisions. Aside from an amplitude change 
Hijing shows little change in correlation structure from 
N-N to central Au-Au collsions. Those results can be 
contrasted with examples from an analysis of real RHIC 
data [1^1 shown in Fig. [S] 




FIG. 5: 2D p t angular autocorrelations from Hijing central 
Au-Au collisions at 200 GeV for quench on (left panel) and 
quench off (right panel) simulations. 

In N-N or p-p collisions there is no obvious quadrupole 
component. However, that possibility is not quanti- 
tatively excluded by the data and requires careful fit- 
ting techniques. For correlation structure which is not 
sinusoidal there is no point to invoking the Wiener- 
Khintchinc theorem on 0a and transforming to a power 
spectrum. Instead, model functions specifically suited to 
such structure (e.g., ID and 2D gaussians) are more ap- 
propriate. The power-spectrum model is appropriate for 
some parts of pa, depending on the collision system. We 
consider hybrid decompositions in Sec. IVIIl 
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D. Comparison with conventional methods 

We can compare the power of the autocorrelation tech- 
nique with the conventional approach based on EP es- 
timation in single-particle n bins. There are concerns 
in the single-particle approach about bias or distortion 
from low bin multiplicities. In contrast, the autocorrela- 
tion method is applicable to any collision system with any 
multiplicity as long as some pairs appear in some events. 
The method is minimally biased, and the statistical error 
in ApA I y/Pref is determined only by the number of bins 
and the total number of events in an ensemble. 

There is interest within the conventional context in 
"correlating" flow results in different n bins. "A three- 
dimensional event shape can be obtained by correlat- 
ing and combining the Fourier coefficients in different 
longitudinal windows" [2l|. That goal is automatically 



achieved with the joint autocorrelation technique but has 
not been implemented with the conventional method. 
The content of an autocorrelation bin at some t?a rep- 
resents a covariance averaged over all bin pairs at n a and 
rjb satisfying t/a = rj a — r\b- If a flow component is iso- 
lated (e.g., V 2 2 /n), a given autocorrelation element rep- 
resents normalized covariance n a n&uf b j 'y/n a nb: where 
vfnab wou ld be the event- wise result of a 'subevent' or 
'scalar-product' analysis between bins a and b. Averaged 
over many events one obtains good resolution on rapid- 
ity and azimuth for arbitrary structures, without model 
dependence or bias. 



VI. FLOW AND MINUETS (NONFLOW) 

The relation between azimuth multipoles and minijets 
is a critically important issue in heavy ion physics which 
deserves precise study. We should carefully compare 
multipole structures conventionally attributed to hydro- 
dynamic flows and parton fragmentation dominated by 
minijets in the same analysis context. The best arena for 
that comparison is the 2D (joint) angular autocorrelation 
and corresponding power-spectrum elements. Before pro- 
ceeding to autocorrelation structure we consider nonflow 
in the conventional flow context 



A. Nonflow and conventional flow analysis 

In conventional flow analysis azimuth correlation struc- 
ture is simply divided into 'flow' and 'nonflow,' where 
the latter is conceived of as non-sinusoidal structure of 
indeterminant origin. The premise is that all sinusoidal 
structure represents flows of hydrodynamic origin. It is 
speculated that nonflow is due to resonances, HBT and 
jets, including minijets. Various properties are assigned 
to nonflow which are said to distinguish it from flow [26| . 
Nonflow is by definition non-sinusoidal and is not corre- 
lated with the RP, thus it can appear perpendicular to 
the RP. 



The multiplicity dependence of nonflow is said to be 
quite different from flow, where "multiplicity depen- 
dence" can sometimes be read as centrality dependence. 

For instance, Q ma ■ Q mb = V rna V mb cos(* ma - \& m6 ) if A, 
B are disjoint, since there are no self pairs. The ensemble 
average then measures covariance V^ ab . For to = 2 it is 
claimed that the nonflow component of V 2 2 ab is cx n c [13] • 
Therefore, V 2 /n cx c, a constant for nonflow — no cen- 
trality dependence. But Vj^/n cx Ap^/ y/p re f which, for 
the minijets dominating nonflow, is very strongly depen- 
dent on centrality [9|, Il2j : the conventional assumption 
is incorrect. The above notation is inadequate because 
minijets (nonflow) should not be included in the Fourier 
power spectrum. They should be modeled by different 
functional forms which wc consider in the next section. 



B. Cumulants 

Another strategy for isolating flow from nonflow is to 
use higher cumulants [2?| . The basic assumption (a phys- 
ical correlation model) is that flow sinusoids are collec- 
tive phenomenon characteristic of almost all particles, 
whereas nonflow is a property only of pairs, termed "clus- 
ters." That scenario is said to imply that v m should be 
the same no matter what the multiplicity, whereas non- 
flow should fall off as some inverse power of n. 

For instance, by subtracting v 2 [A] (four-particle cumu- 
lant) from v 2 [2] (two-particle cumulant) one should ob- 
tain "nonflow" as the difference (cf. Eq. (10) of [Ii|). In 
Fig. 3 1 of we find a plot of g 2 = N part (y\ [2] - v\ [4] ) cx 
Np ar t/n ch x Ap/^/p^J. Multiplying g 2 by n 2 ch /N 2 art we 
obtain a measure of minijet correlations per participant 
pair. That 'nonflow' component increases rapidly with 
centrality (and therefore n), consistent with actual mea- 
surements of minijet centrality trends. The incorrect fac- 
tor in the definition of 172 removes a factor 2x increase 
from peripheral to central in the minijet centrality trend, 
thus suppressing the centrality dependence of 'nonflow.' 



C. Counterarguments 

The conventional flow analysis method requires a com- 
plex strategy to distinguish flow from nonflow in projec- 
tions onto ID azimuth difference </>a- An intricate and 
fragile system results, with multiple constraints and as- 
sumptions. The assumptions are not a priori justified, 
and must be tested. 'Flow' isolated with those assump- 
tions can and does contain substantial systematic errors. 

Claims about the multiplicity (centrality) dependence 
of 'nonflow' (independent of centrality or slowly varying) 
are unsupported speculations without basis in experiment. 
In fact, detailed measurements of minijet centrality de- 
pendence [§, are quite inconsistent with typical as- 
sumptions about nonflow. Multiplicity (centrality) de- 
pendence of flow measurements is further compromised 
by biases resulting from improper statistical methods, 
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especially true for small multiplicities or peripheral colli- 
sions. Such biases can masquerade as physical phenom- 
ena. 

Finally, it is assumed that nonflow has no correlation 
with the RP, thus implying the ability of and need for 
the EP to distinguish flow from nonflow. But nonflow 
(minijets) should be strongly correlated with the EP (jet 
quenching), and such correlations should be measured. 
That is the main subject of paper II in this two-part se- 
ries. Precise decomposition of angular correlations into 
'flow' sinusoids and minijet structure is realized with 2D 
joint angular autocorrelations combined with proper sta- 
tistical techniques. 

Non-flow is a suite of physical phenomena, each wor- 
thy of detailed study. In the conventional approach this 
physics is seen in limited ways by various projections and 
poorly-designed measures and described mainly by spec- 
ulation. With more powerful analysis methods it is pos- 
sibl to separate flow from the various sources of 'nonflow'' 
reliably and identify those sources as interesting physical 
phenomena. 



VII. STRUCTURE OF THE JOINT ANGULAR 
AUTOCORRELATION IN A-A COLLISIONS 

We now return to the 2D angular autocorrelation. By 
separating its structure into a few well-defined compo- 
nents we obtain an accurate separation of multipoles, 
minijets and other phenomena. Minijets and "flows" can 
be compared quantitatively within the same analysis con- 
text. 

Each bin of an autocorrelation is a comparison of two 
"subevents." The notional term "subevent" represents a 
partition element in conventional math terminology (e.g., 
topology, cf. Borel measure theory). An "event" is a dis- 
tribution in a bounded region of momentum space (de- 
tector acceptance), and a subevent is a partition element 
thereof. A distribution can be partitioned in many ways: 
by random selection, by binning the momentum space, 
by particle type, etc. A uniform partition is a binning, 
and the set of bin entries is a histogram. 

A bin in an angular autocorrelation represents an av- 
erage over all bin pairs in single-particle space separated 
by certain angular differences (t)a, 0a)- The bin contents 
represent normalized covariances averaged over all such 
pairs of bins. The notional "scalar product method," re- 
lating two subevents in conventional flow analysis, is al- 
ready incorporated in conventional mathematical meth- 
ods developed over the past century as covariances in 
bins of an angular autocorrelation. Using 2D angular au- 
tocorrelations we easily and accurately separate nonflow 
from flow. "Nonflow" so isolated has revealed the physics 
of minijets — hadron fragments from the low-momentum 
partons which dominate RHIC collisions. 



A. Minijet angular correlations 

Minijet correlations are equal partners with multi- 
pole correlations on difference- variable space (t]a,4>a.)- 
Minimum-bias jet angular correlations (dominated by 
minijets) have been studied extensively for p-p and Au- 
Au collisions at 130 and 200 GeV [! EES El • Those 
structures dominate the "nonflow" of conventional flow 
analysis. In p-p collisions minijet structure — a same-side 
peak (jet cone) and away-side ridge uniform on rj^ — 
are evident for hadron pairs down to 0.35 GeV/c for 
each hadron. Parton fragmentation down to such low 
hadron momenta is fully consistent with fragmentation 
studies over a broad range of parton energies (e.g., LEP, 
HERA) 




FIG. 6: 2D pt angular autocorrelations from Au-Au collisions 
at 200 GeV for 80-90% central collisions (left panels) and 45- 
55% central collisions (right panels). In the lower panels si- 
nusoids cos(</>a) and cos(2<^a) have been subtracted to reveal 
"nonflow" structure. 

In Fig. [6] we show autocorrelations obtained by inver- 
sion of pt fluctuation scale (bin-size) dependence [9( . The 
upper-left panel is 80-90% central and the upper-right 
panel is 45-55% central Au-Au collisions. Correlation 
structure is dominated by a same-side peak and mul- 
tipole structures (sinusoids). Subtracting the sinusoids 
reveals the minijet structure in the bottom panels and 
illustrates the precision with which flow and nonflow can 
be distinguished. The negative structure surrounding the 
same-side peak at lower right is an interesting and unan- 
ticipated new feature Q. 

B. Decomposing 2D angular autocorrelations: a 
controlled comparison 

Based on extensive analysis [H, we find 

three main contributions to angular correlations in RHIC 
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nuclear collisions: 1) transverse fragmentation (mainly 
minijets), 2) longitudinal fragmentation (modeled as 
"string" fragmentation), 3) azimuth multipoles (flows). 
Longitudinal fragmentation plays a reduced role in heavy 
ion collisions. In this study we focus on the interplay be- 
tween 1) and 3), transverse parton fragmentation and 
azimuth multipoles, as the critical analysis issue for az- 
imuth correlations in A-A collisions. 

The 2D joint autocorrelation Pa(va,^a) is the basis 
for decomposition. The criteria for distinguishing az- 
imuth multipoles from minijet structure are n& depen- 
dence and sinusoidal 4>a dependence. Structure with si- 
nusoidal 4>a dependence and tja invariance is assigned to 
azimuth multipoles. Other structure, varying generally 
on (ijaj^a), is assigned in this exercise to minijets. We 
adopt the decomposition 



Pa(va,4>a) = Pj(r)A,4>A) +Pm(<M), 



(55) 



where j represents (mini)jets and m represents multi- 
poles. That decomposition is reasonable within a limited 
pseudorapidity acceptance, e.g., the STAR TPC accep- 
tance [29|. Over a larger acceptance other separation 
criteria must be added. 

To illustrate the separation process we construct artifi- 
cial autocorrelations combining flow sinusoids and mini- 
jet structure with centrality dependence taken from mea- 
surements. We add statistical noise appropriate to a typ- 
ical event ensemble of a few million Au-Au collisions. We 
then fit the autocorrelations with model functions and x 2 



minimization. We compare the resulting fit parameters 
with the input parameters. We then project the 2D au- 
tocorrelations onto 4>a and fit the results with a sinusoid. 
The result of the ID fit represents the product of a con- 
ventional flow analysis if the EP resolution is perfectly 
corrected and there is no statistical bias in the method. 
We then compare the resulting sinusoid amplitudes. 

In Fig. [7]we show an analysis for 80-90% central Au-Au 
collisions, which are nearly N-N collisions. The upper- 
right panel is an accurate model of data from p-p col- 
lisions 11]. The upper-left panel is the simpler repre- 
sentation for this exercise with added statistical noise. 
The difference in constant offsets is not relevant to the 
exercise. At lower left is the distribution "seen" by a 
conventional flow analysis, which simply integrates over 
(projects) the r\ dependence. The projection includes a 
triangular acceptance factor on r)& imposed on the joint 
autocorrelation, resulting in distortions that affect the si- 
nusoid fit. The lower-right panel is the projection onto 
<f) A with x 2 fit. The ID fit gives Ap[2]/^p^J = 0.0052, 
compared to 0.0013 from the 2D fit 

In Fig. [8] we show the same analysis applied to 5-10% 
Au-Au collisions. The model distribution, derived from 
observed data trends, is dominated by a dipole term 
oc cos(^a) and an elongated same-side jet peak, although 
the combination closely mimics a quadrupole oc cos(20a) 
(elliptic flow) . The lower- left panel shows the effect of the 
ry-acceptance triangle on tja applied to the upper-right 
panel which is implicit in any projection onto 4>a (obvious 
in this 2D plot). The resulting projection on </>a is shown 





FIG. 7: Simulated three-component 2D angular autocorrela- 
tion for 80-90% central Au-Au collisions at 200 GeV (upper 
left), model of data distribution from fitting (upper- right), 
autocorrelation with eta-acceptance triangle imposed (lower 
left) and sinusoid fit to ID projection of lower-left panel 
(lower-right). 




« 0.4 F 




FIG. 8: Same as the previous figure but for 5-10% central Au- 
Au collisions at 200 GeV. Note the pronounced effect of the 
^-acceptance triangle in the lower-left panel (relative to the 
upper-right panel) resulting from projection of space (^1,772) 
onto its difference axis. 
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at lower right. The ID fit gives Ap[2]/^/p7^7 = 0.103, 
compared to 0.060 from the 2D fit. The differences be- 
tween ID and 2D fits are much larger than the differences 
between input and fitted parameters in the previous ex- 
ercise. They reveal the limitations of conventional flow 
analysis. 

In Fig. we give a summary of results for twelve cen- 
trality classes, nine 10% bins with mean values 95% • • • 
15%, plus two 5% bins with mean values 7.5% and 2.5%. 
The twelfth class is b — (0%), constructed by extrapo- 
lating the parameterizations from data. The solid curves 
and points represent the parameters inferred from 2D fits 
to the joint angular autocorrelations. The dashed curves 
represent the model parameters used to construct the 
simulated distributions. There is good agreement at the 
percent level. 




FIG. 9: A parameter summary for the previous two figures. 
The solid curves are input model parameters for m = 1 and 
m = 2 sinusoids cos(m</>A) and a same-side 2D gaussian with 
two widths and peak amplitude which approximate 200 GeV 
Au-Au collisions. The dashed curves are the results of fits to 
the model 2D autocorrelations exhibiting excellent accuracy. 
The dash-dot curve represents ID fits to projections on <pA 
(previous lower-right panels) corresponding to conventional 
flow analysis. 

The fits to ID projections on 4>a (dash-dot curve) how- 
ever differ markedly from the 2D fit results and the in- 
put parameters. The differences are very similar to the 
changes of conventional flow measurements with different 
strategies to eliminate "nonflow." This exercise demon- 
strates that with the 2D autocorrelation there is no guess- 
work. We can distinguish the multipole contributions 
from the minijet contributions. The 2D angular autocor- 
relation provides precise control of the separation. 



C. V2 in various contexts 

In Fig. [10] we contrast the results of ID conven- 
tional flow analysis (dashed curves) and extraction of the 
quadrupole amplitude from the 2D angular autocorrela- 
tion (solid curves). We make the correspondence 



among variables, with n/2i: modeling po = cPn/drjdcj). 
That expression defines a minimally-biased vi . We make 
the comparison in four plotting formats. The upper-right 
panel is most familiar from conventional flow analysis, 
where v-i is plotted vs participant nucleon number. The 
trends can be compared with Fig. 13 of [3(| (open circles 
vs solid stars). Also included in that panel is the trend 
V2 ~ 0.22 e predicted by hydro for thermalized A- A col- 
lisions at 200 GeV (dotted curve, [U). 
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FIG. 10: The quadrupole component in various plotting con- 
texts derived from the previous model exercise. Conventional 
flow measure V2 is shown in the upper panels. Ap[2]/ y/pref 
based on Pearson's normalized covariance is shown in the 
lower panels. The upper-right panel shows «2 vs n par t in the 
conventional plotting format. The lower-left panel repeats 
the left panel of the previous figure, with conventional ID fit 
results shown as the dashed curve. Dashed and solid curves 
correspond in all panels, v estimates the mean N-N encoun- 
ters per participant pair. The dotted curve at lower right is 
the error function used to generate the model quadrupole am- 
plitudes. It's correspondent for V2 is at upper left. The dotted 
curve at upper right is eccentricity e from a parameterization. 

The remaining panels are plotted on parameter v = 
Znbin/npart, the ratio of N-N binary encounters to par- 
ticipant pairs which estimates the mean participant path 
length in number of encountered nucleons, a geometri- 
cal measure. Comparing the upper panels we see that v 
treats peripheral and central collisions equitably, whereas 
npart or n c harge compresses the important peripheral re- 
gion into a small interval. 

In the lower-left panel we plot per-particle density ra- 
tio A /cm [2]/ y/Pref vs v. That quantity, when extracted 
from 2D fits, rises from near zero for peripheral collisions 
to a maximum for mid-central collisions, falling toward 
zero again for 6 = 0. In contrast to V2 , which is the 
square root of a per-pair correlation measure, the per- 
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particle density ratio reflects the trend for "flow" in the 
sense of a current density. "Flow" is small for periph- 
eral collisions and grows rapidly with increasing nucleon 
path length. The trend with centrality is intuitive. The 
values obtained from the ID projection per conventional 
flow analysis (dashed curve) are consistently high, espe- 
cially for central collisions, exhibiting a strong systematic 
bias. The ID fit procedure (identical to the "standard" 
and two-particle flow methods) confuses minijet structure 
with quadrupole structure. 

In the lower-right panel we show the density ratio di- 
vided by initial spatial eccentricity e defined by a pa- 
rameterization derived from a Glauber simulation and 
plotted as the dotted curve in the upper-right panel [3l[ . 
The trend from 2D fits (solid curve) is closely approx- 
imated by a simple error function (dotted curve) with 
half-maximum point at the center of the v range. In 
fact, the dotted curve is the basis for generating the input 
quadrupole amplitudes for our model, and the small devi- 
ation of the solid curve from the dotted curves in upper- 
left and lower-right panels reveals the systematic error or 
bias in the 2D fitting procedure (~ 20% at v ~ 1, < 5% 
at v = 5.8). The dashed curves from the conventional ID 
fits show large relative deviations from the input trend, 
especially for peripheral collisions where the emergence of 
collectivity is of interest and for central collisions where 
the issue of thermalization is most important. 

Comparing the solid curve to existing data in the 
format of the upper-right panel shown in [301 ] (Fig. 13) 
indicates that our simple formulation at lower right 
(dotted curve) is roughly consistent with analysis of 
flow data based on four-particle cumulants. The re- 
sult from data generated with the same model and ana- 
lyzed with the conventional flow analysis method (dashed 
curve in upper-right panel) also agrees with the conven- 
tional method applied to real RHIC data. Our model 
may therefore indicate an underlying simplicity to the 
quadrupole mechanism which is not hydrodynamic in ori- 
gin. The model centrality trend for v-i is certainly incon- 
sistent with the hydro expectation vi cx e (20l . [23|, as 
demonstrated in the upper-right panel (c/. App. [E)) . 



VIII. DISCUSSION 

A. Conventional flow analysis 

The overarching premise of conventional flow analy- 
sis is that in the azimuth distribution of each collision 
event lies evidence of collective phenomena which must 
be discovered to establish event-wise thermalization. The 
1/n — *■ hydro limit shapes the analysis strategy, and 
flow manifestations are the principal goal. Finite event 
multiplicities are seen as a major source of systematic 
error, as are correlation structures other than flow. Mul- 
tiple strategies are constructed to deal with non-flow and 
finite multiplicities. A stated advantage of the conven- 
tional method is that the Fourier coefficients can be cor- 



rected. The great disadvantage is they mustbe corrected. 

From the perspective of two-particle correlation anal- 
ysis, especially in the context of autocorrelations and 
power spectra, the conventional program leaves much to 
be desired. The conventional analysis is essentially an 
attempt to measure two-particle correlations with single- 
particle methods combined with RP estimation, similar 
to the use of trigger particles in high-;^ jet analysis. But, 
by analysis of the algebraic structure of conventional flow 
analysis we have demonstrated that RP estimation does 
not matter to the end result. Without a proper statistical 
reference conventional analysis results contain extrane- 
ous contributions from the statistical reference which arc 
partially 'corrected' in a number of ways. The improper 
treatment of random variables incorporates sources of 
multiplicity-dependent bias in measurements, and the fi- 
nal results are questionable. 

Flow measure v m is nominally the square root of per- 
pair correlation measure V^/n(n— 1). v m centrality 
trends are thus nonintuitive and misleading (e.g., "el- 
liptic flow" decreases with increasing A-A centrality). 
The situation is similar to per-pair fluctuation measure 
S pt , which provides a dramatically misleading picture of 
Pt fluctuation dependence on collision energy In 
contrast, per-particle correlation measures provide intu- 
itively clear results and often make dynamical correla- 
tion mechanisms immediately obvious. In particular, the 
mechanisms behind "nonflow" in the form of minijets are 
clearly apparent when correlations are measured by per- 
particle normalized covariance density Ap/ ^p re f in a 2D 
autocorrelation. 



B. Autocorrelations and nonflow 

In conventional flow analysis it is proposed to measure 
flow in narrow rapidity bins ("strips") so as to develop 
a three-dimensional picture of event structure. However, 
there has been little implementation of that proposal. In 
contrast, the joint angular autocorrelation by construc- 
tion contains all possible covariances among pseudora- 
pidity bins within a detector acceptance. The ideal of 
full event characterization is thereby realized. 

The angular autocorrelation is the optimum solution to 
a geometry problem — how to reduce the six-dimensional 
two-particle momentum space to two-dimensional sub- 
spaces with minimum distortion or information loss. The 
autocorrelation is the unique solution to that problem, 
involving no model assumptions. The ensemble-averaged 
angular autocorrelation contains all the correlation infor- 
mation obtainable from a conventional flow analysis, but 
with negligible bias and no sensitivity to individual event 
multiplicities. 

Because it is a two-dimensional representation the 
angular autocorrelation is far superior for separating 
"flow" (multipoles) from "nonflow" (minijets), as we have 
demonstrated. A simple exercise demonstrates that sep- 
aration is complete at the percent level, whereas the con- 
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ventional method admits crosstalk at the tens of percent 
level. Precise separation leads to new physics insights 
from the multipoles and minijets so revealed. 



C. Collision centrality dependence 

Collision centrality dependence is of critical impor- 
tance in the comparison of flow and minijets. Parton col- 
lisions and hydrodynamic response to early pressure have 
very different dependence on impact parameter and col- 
lision geometry, especially for peripheral collisions. Pe- 
ripheral A-A collisions should approach p-p (N-N) col- 
lisions, and correlation structure may change rapidly in 
mid-peripheral collisions as collective phenomena develop 
there. The possible onset of collective behavior in mid- 
peripheral collisions and reduction in more central col- 
lisions are of major importance for understanding the 
relation of minijets to flow. The conventional flow analy- 
sis method is severely limited for peripheral collisions. 
In contrast, correlation measure A.p/^/p re f, centrality 
measure v and associated centrality techniques described 
in [321 ] are uniquely adapted to cover all centrality regions 
down to N-N with excellent accuracy. 



D. Physical interpretations 

Because similar flow measurement techniques have 
been applied at Bevalac and RHIC energies with sim- 
ilar motivations it is commonly assumed that azimuth 
multipoles have a common source over a broad collision 
energy range — hydrodynamic flows, collective response 
to early pressure. The hydro mechanism was proposed 
as the common element in [20[ and persists as the lone 
interpretation of azimuth multipoles in HI collisions to 
date. 

At Bevalac and AGS energies it is indeed likely that 
azimuth multipoles result from 'flow' of initial-state nu- 
cleons in response to early pressure, with consequent 
final-state correlations of those nucleons — a true hydro 
phenomenon. However, at SPS and RHIC energies the 
source of azimuth multipoles inferred from final-state 
produced hadrons (mainly pions) may not be hydrody- 
namic, in contrast to arguments by analogy with lower 
energies. Other sources of multipole structure should be 
considered [34], [35| • Multipoles at higher energies could 
arise at the partonic or hadronic level, early or late in the 
collision, with collective motion or not, and if collective 
then implying thermalization or not. 

The chain of argument most often associated with 
elliptic flow asserts that observation of flow as a col- 
lective phenomenon demonstrates that a thermalized 
medium (QGP) has been formed which responds hy- 
drodynamically to early pressure and converts an ini- 
tial configuration-space eccentricity to a corresponding 
quadrupole moment in momentum space. 



However, nonflow in the form of minijets provides con- 
tradictory evidence. Minijet centrality trends indicate 
that thermalization is incomplete, and substantial mani- 
festations of initial-state parton scattering remain at ki- 
netic decoupling 0, [l(| EH- Precision studies of mini- 
jet centrality dependence (y dependence) indicate that a 
large fraction of the minijet structure expected from lin- 
ear superposition of N-N collisions (no thermalization) 
persists in central Au-Au collisions. That contradic- 
tion requires more complete experimental characteriza- 
tion and careful theoretical study (3o| . 

Arguments based on interpreting the quadrupole com- 
ponent as hydrodynamic flow exclude alternative phys- 
ical mechanisms. Aside from minijet systematics there 
are other hints that a different mechanism might be re- 
sponsible for azimuth multipoles. In Fig. [10] we showed 
that flow measurements based on four-particle cumulants 
(with bias sources and nonflow thereby reduced) are best 
described by a trend (solid curves) that is inconsistent 
with the hydro expectation V2 oc e. The trend is in- 
stead simply described in terms of per-particle measure 
Ap/^/p re f and two shape parameters relative to e. 

We question the theoretical assumption that e should 
be simply related to v-i as opposed to some other mea- 
sure of the azimuth quadrupole component. We expect 
a priori and find experimentally that variance measures, 
integrals over two-particle momentum space, more typ- 
ically scale linearly with geometry parameters. Thus, 
Ap[2]/ V ^^y cx nv\ may be more closely related to e, 
and the relation may or may not be characteristic of a 
hydro scenario. 



IX. SUMMARY 

In conclusion, we have reviewed Fourier transform the- 
ory, especially the relation of autocorrelations to power 
spectra, essential for analysis of angular correlations in 
nuclear collisions. In that context we have reviewed 
five papers representative of conventional flow analysis 
and have related the methods and results to autocorre- 
lation structure and spherical and cylindrical multipole 
moments. 

We have examined the need for event-plane evaluation 
in correlation measurements and find that it is extraneous 
to measurement of azimuth multipole moments. The EP 
estimate drops out of the final ensemble average. 

We have introduced the definition of the 2D (joint) an- 
gular autocorrelation and considered the distinction be- 
tween flows (cylindrical multipoles) and nonflow (domi- 
nated by minijet structure) in conventional flow analysis 
and criticized the basic assumptions used to distinguish 
the two in that context. 

Based on measured minijet and flow centrality trends 
we have constructed a simulation exercise in which model 
autocorrelations of known composition are combined 
with statistical noise from a typical event ensemble and 
fit with a model function consisting of a few simple com- 
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ponents, first as a 2D autocorrelation and second as a ID 
projection on azimuth difference axis 4>\. We show that 
the 2D fit returns input parameters accurately at the 
percent level, whereas the ID fit, representing conven- 
tional flow analysis, deviates systematically and strongly 
from the input. Comparisons with published flow data 
indicate that the observed bias in the simulation is ex- 
actly the difference attributed to "nonflow" in conven- 
tional measurements. 

By comparing our simple algebraic model of 
quadrupole centrality dependence to data we observe 
that the trend vi oc e is not met for any collision sys- 
tem, nor is there asymptotic approach to such a trend. 
That observation raises questions about the relevance of 
hydrodynamics to phenomena currently attributed to el- 
liptic flow at the SPS and RHIC. 

This work was supported in part by the Office of Sci- 
ence of the U.S. DoE under grant DE-FG03-97ER41020. 



APPENDIX A: BROWNIAN MOTION 

There is a close analogy between Brownian motion and 
the azimuth structure of nuclear collisions. The long his- 
tory of Brownian motion and its mathematical descrip- 
tion can thus provide critical guidance for the analysis 
of particle distributions. Brownian motion (more gen- 
erally, random motion of particles suspended in a fluid) 
was modeled by Einstein as a diffusion process (random 
walk) [7]. He sought to test the "kinetic- molecular" the- 
ory of thermodynamics and provide direct observation of 
molecules. Paul Langevin developed a differential equa- 
tion to describe such motion, which included a stochastic 
term representing random impulses delivered to the sus- 
pended particle by molecular collisions. Jean Perrin and 
collaborators performed extensive measurements which 
confirmed Einstein's predictions andprovided definitive 
evidence for the reality of molecules Q. 



quasi-random processes. Measuring multipole moments 
on azimuth in nuclear collisions is formally equivalent 
to measuring "drift" terms on time in the quasi-random 
walk of a charged particle suspended in a molecular fluid 
within a superposition of oscillating electric fields. There 
are many other applications for Eq. (|A1|) . 

For a true random walk consisting of uncorrelated steps 
Einstein expressed (r 2 ) /8t = d-2D (random walk in d di- 
mensions) in terms of diffusion coefficient D. The second 
term (r 2 cos(</>a)) = (St) 2 v 2 represents a possible deter- 
ministic component (correlations) , with x the direction of 
an applied "force." In that case successive angles <fii are 
correlated, and the result is a macroscopic nonstochastic 
drift of the walker trajectory. 

The fractal dimension of a random walk [first term 
in Eq. (|A1[) ] is df = 2. The trajectory is therefore 
a "space-filling" curve in 2D configuration space. The 
appropriate measure of trajectory size is area, and the 
rate of size increase is the diffusion coefficient (rate of 
area increase). In contrast, the second term in Eq. (|A1|) 
represents a deterministic trajectory whose nominal di- 
mension is 1 (modulo the extent of curvature, which in- 
creases the dimension above 1). Therefore, the appropri- 
ate measure of trajectory size is length, and speed is the 
correct rate measure. For Brownian motion with drift 
the trajectory dimension is not well-defined, depending 
on the relative magnitudes of the drift and stochastic 
terms, and the concept of speed is therefore ambiguous. 
Attempts to measure the linear speed of Brownian mo- 
tion in the nineteenth century failed because of the frac- 
tal structure of random walks. From the structure of 
Eq. (|A1|) the average speed over interval At — ndt is 
^/R 2 /(At) 2 ~ (r 2 /{St) 2 )/n, and the limiting case for 
At = nSt — > is the so-called "infinite speed of diffu- 
sion." That topological oddity is formally equivalent to 
the "multiplicity bias" of conventional flow analysis. 

2. Brownian motion and nuclear collisions 



1. The quasi-random walker 

We model a 2D quasi-random walker (including 
nonzero correlations) as follows. The walker position 
is recorded in equal time intervals St. After n steps, 
with step-wise displacements r sampled randomly from a 
bounded distribution, the walker position relative to an 
arbitrary starting point is, in the notation of this paper, 
R = where r% is the i th displacement. The 

squared total displacement is then 

R 2 = n(r 2 ) +n(n- l)(r 2 cos(</> A ))- (Al) 

The first term, linear in n (or t), was described by Ein- 
stein. The second term could represent "drift" of the 
walker due to deterministic response to an external in- 
fluence. The composite is then termed "Brownian motion 
with drift," a popular model for stock markets and other 



We now consider the close analogy between Einstein's 
theory of Brownian motion and the measurement of p 2 
in a nuclear collision, using directivity as an example. 
Just as R is the vector total displacement of a quasi- 
random walker in 2D configuration space, Q\ is the vector 
total displacement of a quasi-random walker (event-wise 
particle ensemble) in 2D momentum space. After n steps 
the squared displacements are 

R 2 = n 2 St 2 v' 2 = ndt AD + n(n - 1) 5t 2 v 2 x (A2) 
Ql = n 2 p' 2 = n{p 2 )+n(n-l)p 2 x . 

AD5t is the increase in area per step of a random walker 
in 2D configuration space, (p 2 ) is the increase in area per 
step (per particle) of a random walker in 2D momentum 
space, playing the same role as the diffusion coefficient. 
The RHS first term in the first line is the subject of Ein- 
stein's 1905 Brownian motion paper. Its measurement by 
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Perriii confirmed the reality of molecules and the validity 
of Boltzmann's kinetic theory. 

As noted, attempts to measure mean speed v' x of a par- 
ticle in a fluid failed because speed is the wrong rate mea- 
sure for trajectory size increase. Speed measurements 
decreased with increasing sample number or observation 
time. It was not until Einstein's formulation and later 
mathematical developments that the topology of the ran- 
dom walk and its consequences became apparent. Initial 
attempts at the Bevalac to measure p x in the form p' x 
using directivity failed for the same reason. Corrections 
were developed to approximate the unbiased quantity p x , 
and the failure was attributed to multiplicity bias or 'au- 
tocorrelations.' Ironically, the autocorrelation distribu- 
tion is the ideal method to access the unbiased quantity 
in either case. 



3. Einstein and autocorrelations 

To provide a statistical description of Brownian motion 
Einstein introduced the autocorrelation concept with the 
following language Q- 

Another important consideration can be re- 
lated to this method of development. We 
have assumed that the single particles are all 
referred to the same co-ordinate system. But 
this is unnecessary, since the movements of 
the single particles are mutually independent. 
We will now refer the motion of each parti- 
cle to a co-ordinate system whose origin co- 
incides at the [arbitrary] time t = with the 
[arbitrary] position of the center of gravity of 
the particle in question; with this difference, 
that [probability distribution] f(x,t)dx now 
gives the number of the particles whose x co- 
ordinate has increased between the time t = 
and the time t = t, by a quantity which lies 
between x and x + dx. 

Einstein's function /(£, t) is a 2D autocorrelation 
which satisfies the diffusion equation. The solution is 
a gaussian on x relative to an arbitrary starting point 
(thus defining difference variables £ = x — x s t a rt and 
t = t — t s tart), with ID variance <r| = 2Dt. The au- 
tocorrelation is sometimes called a two-point correlation 
function or two-point autocorrelation. The angular auto- 
correlation is a wide-spread and important analysis tool, 
e.g., in astrophysics, nuclear collisions and many other 
fields. 



4. Wiener, Khintchine, Levy and Kolmogorov 

The names Wiener, Levy, Kolmogorov and Khintchine 
figure prominently in the copious mathematics derived 
from the Brownian motion problem. Norbert Wiener led 



efforts to provide a mathematical description of Brown- 
ian motion, abstracted to a Wiener process, a special case 
of a Levy process (generalization of a discrete random 
walk to a continuous random process) [37j • The Wiener- 
Khintchinc theorem provides a power-spectrum represen- 
tation for stationary stochastic processes such as random 
walks, for which a Fourier transform does not exist. We 
have acknowledged the theorem with our Eq. (|10[) . 

The analysis of azimuth structure in nuclear collisions 
in terms of angular autocorrelations is based on power- 
ful mathematics developed throughout the past century. 
Autocorrelations make it possible to study azimuth struc- 
ture for any event multiplicity down to p-p collisions with 
as little as two detected particles per event. The effects of 
"non-flow" can be eliminated from "flow" measurements 
(and vice versa) without model dependence or guesswork. 
The Brownian motion problem and Einstein's fertile so- 
lution inform two central issues for studies of the correla- 
tion structure of nuclear collisions: analysis methodology 
and physics interpretation. 



APPENDIX B: RANDOM VARIABLES 



A random variable represents a set of samples from 
a parent distribution. The outcome of any one sample 
is unpredictable (i.e., random), but through statistical 
analysis an ensemble of samples can be used to infer prop- 
erties (statistics - results of algorithms applied to a set of 
samples) of the parent distribution. Sums over particles 
and particle pairs of kinematic quantities are the primary 
random variables in analysis of nuclear collision data. 



1. The algebra of random variables 

Products and ratios of random variables behave non- 
intuitively because random variables don't obey the alge- 
bra of ordinary variables. E.g., factorization of random 
variables results in the spawning of covariances. The 
approximation xy ~ xy common in conventional flow 
analysis is a source of systematic error (bias) because 
xy = xy + xy — xy. The omitted term is a covariance. 
Such covariances play a role in statistics similar to QM 
commutators, with 1/n <-» K. Conventional flow analysis 
assumes the 1/n — > limit for some random variables, 
and the results are undependable for small multiplicities. 
Similarly, improper treatment of ratios of random vari- 
ables results in infinite series of covariances. E.g., 



— — x , ox ■ on x ■ yon)" 
x/n=-(l+ + + — h^- + •••), Bl 

n xn xn^ 

with (6n) 2 /n = cr^/n ~ 1 — 2. Thus, the common ap- 
proximation x/n ~ xjh can result in significant n- and 
physics-dependent (x-n covariances) bias for small n. 



22 



In this paper we distinguish between event-wise and 
ensemble-averaged quantities and do not employ en- 
semble averages of ratios of random variables. We in- 
clude event-wise factorizations and ratios only to sug- 
gest qualitative connections with conventional flow anal- 
ysis. E.g., we consider V£ = n{n — 1)(cos(to</>a)) with 
(cos (m</> A )) = (cos 2 (m[^ - * r ])) = v 2 m . But, ^ 
V^jn(n — 1) 7^ vf n . v m as typically invoked in conven- 
tional flow analysis is not a well-defined statistic. 

2. Statistical references 

The concept of a statistical reference is largely absent 
from conventional flow analysis. By 'statistical reference' 
we mean a quantity or distribution which represents an 
uncorrelated system, a system consistent with indepen- 
dent samples from a fixed parent distribution (central 
limit conditio ns l al). Concerns about 'bias' from low mul- 
tiplicities [III [21U22II typically relate to the presence of an 
unsubtracted and unacknowledged statistical reference 
in the final result. Finite multiplicity fluctuations are 
then said to produce systematic errors, false azimuthal 
anisotropics, a problem masking true collective effects. 

In the limit l/n — > the statistical reference may in- 
deed become negligible compared to the true correlation 
structure. However, its presence for nonzero l/n is a po- 
tential source of systematic error which may block access 
to important small-multiplicity systems (peripheral col- 
lisions and/or small kinematic bins). In general, if the 
statistical reference is not correctly subtracted the result 
is increasingly biased with smaller multiplicities. Identi- 
fication and subtraction of the proper reference is one of 
the most important tasks in statistical analysis. 

Use of the term 'statistical' to mean 'uncorrelated' is 
misleading (e.g., 'statistical' vs 'dynamical'). All ran- 
dom variables and their fluctuations about the mean are 
'statistical.' Some random variables and their statistics 
are reference quantities, representing systems that are by 
construction uncorrelated (independent sampling from a 
fixed parent). We therefore label statistical reference 
quantities 'ref,' not 'stat.' 

3. Random variables and Fourier analysis 

In the context of Fourier analysis the basic finite- 
number (Poisson) statistical reference is manifested as 
the delta-function component in the autocorrelation den- 
sity Eq. (|13p and the white- noise constant term n(r 2 ) in 
the event-wise power spectrum. Other reference compo- 
nents may arise from two-particle correlations which are 
not of interest to the analysis (e.g., detector effects) and 
which may be revealed in mixed-pair distributions. A 
clear distinction should always be maintained between 
the reference and the sought-after correlation signal. 

Careful attention to random-variable algebra is es- 
pecially important in a Fourier analysis. The power- 



spectrum elements and autocorrelation density must sat- 
isfy the transform equations both for each event and after 
ensemble averaging. In conventional flow analysis that 
condition is often not satisfied. For instance, and 
satisfy the FT transforms and Wiener-Khintchine the- 
orem before and after ensemble averaging respectively, 
whereas the v m do not. 



4. Minimally-biased random variables 

It is frequently stated in the conventional flow liter- 
ature that flow analysis must insure sufficiently large 
multiplicities. The operating assumption in the design 
of conventional flow methods is the continuum limit 
l/n — > 0, with inevitable bias for smaller multiplicities. 
However, careful reference design and algebraic manipu- 
lation of random variables makes possible precise treat- 
ment of event-wise multiplicities down to n = 1. Some 
statistical measures perform consistently no matter what 
the sample number. The full multiplicity range is essen- 
tial to measure azimuth multipole evolution with central- 
ity down to N-N and p-p collisions, so that A- A "flow" 
phenomena may be connected to phenomena observed in 
elementary collisions and understood in a QCD context. 

Since multiplicity necessarily varies strongly with cen- 
trality, multiplicity-dependent bias in flow measurements 
is unacceptable, and every means should be used to in- 
sure minimally-biased statistics. To achieve that end 
analysis methods must carefully transition from safe 
event-wise factorizations (as featured in this paper) to 
ensemble averages minimally biased for all n. Linear 
combinations of powers of random variables, e.g., vari- 
ances and covariances, satisfy a linear algebra. Such in- 
tegrals of two-particle momentum space are nominally 
free of bias. 



APPENDIX C: MULTIPOLES AND SPHERICITY 

The ID Fourier transform on azimuth is part of a larger 
representation of angular structure. The encompassing 
context is a 2D multipole decomposition on (#, </>) repre- 
sented by the sphericity tensor, with the spherical har- 
monics Y™ as elements. In limiting cases submatrices 
of the sphericity tensor reduce to "cylindrical harmon- 
ics" cos(m</>), part of the ID Fourier representation on 
azimuth. 

The central premise of a multipole representation 
is that the final-state particle angular distribution on 
[9(y z ),(j>] is efficiently represented by a few low-order 
spherical harmonics (SH) Y [ m (6, (f>). At the Bevalac, 
sphericity tensor S containing spherical harmonics Y™ as 
elements was introduced. Directivity Q\, simply related 
to Y%, was employed to represent a rotated quadrupole 
as a dipole pair antisymmetric about the collision mid- 
point. At lower energies (Bevalac, AGS) the quadrupole 
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principal axis may be rotated to a large angle with re- 
spect to the collision axis and Y^ dominates. At higher- 
energies and near midrapidity (0 ~ tt/2) the dominant 
SH is F 2 2 . 

1. Spherical harmonics 

The spherical harmonics are defined as 

r^-^-^rr^D^"*. (CD 

where P ; m (0) is an associated Legendre function [38| . An 
event-wise density on the unit sphere can be expanded 
as 

p(Q) = ^QimY^Sl) (C2) 
I in 




dnY; m (n)p(n) 



i 

= n(Y; m (n)), 

where fi -> {9,(j>) and dft = dcos(9)d<p. The FTs on 
form a special case of those relations when p(f2) is peaked 
near 9 ~ tt/2. The Y; m are orthonormal and complete: 

J dSlY lm {6,<t>)Y Vm ,{e,cf>) = 5 w S mm , (C3) 
^F iro (0, 0)1^(0', <//) = 5(Q-ST). 

2. Multipoles 

The spherical harmonics are model functions for single- 
particle densities on (0, </>). The coefficients of the mul- 
tipole expansion of a distribution are complex spherical 
multipole moments describing 2 l poles and defined as en- 
semble averages of the spherical harmonics over the unit 
sphere weighted by an angular density. 

The following relation is defined by analogy with the 
expansion of an electric potential in spherical harmonics, 
in this case on momentum space p rather than configu- 
ration space f [38j 

J dp WP\ ~ ^2TTT qim ^^- (C4) 

Ira 

The coefficients are the event-wise spherical multipole 
moments 

vim = j P 2 dpdn P l p{p^)YC m {n) (C5) 

n 

i 

= n{p l Y* n (n)). 



Eq. (|C2p is the special case for p restricted to unity (i.e., 
distribution on the unit sphere). 

In general, ?RY™ cx sin m (0) cos(m(/>), and the cos(m<j>) 
are by analogy "cylindrical harmonics" [42j ■ The ensem- 
ble average of a cylindrical harmonic over the unit circle 
weighted by ID density p{4>) results in complex cylindri- 
cal multipole moments Q m . The Fourier coefficients Q m 
obtained from analysis of SPS and RHIC data are there- 
fore cylindrical multipole moments describing 2m poles. 
E.g., m — 2 denotes a quadrupole moment and m = 4 
denotes an octupole moment,. 

If nonflow contributions (i.e., structure rapidly vary- 
ing on r\ or y) are present, a multipole decomposition of 
p(9, <f) is no longer efficient, and the inferred multipole 
moments are difficult to interpret physically (e.g., flow in- 
ferences per se are biased). In Sec. [V] we describe a more 
differential method for representing angular structure us- 
ing two-particle joint angular autocorrelations on differ- 
ence axes (?7A)<M)- Given a decomposition of p{0,4>) 
based on variations on tja we can distinguish cylindri- 
cal multipoles accurately from "nonflow" structure (c/. 
Sec. 

3. Sphericity 

The sphericity tensor has been employed in both jet 
physics and flow studies. A normalized 3D sphericity 
tensor was defined in [14| to search for initial evidence 
of jets in e + -e~ collisions. A decade later sphericity was 
introduced to the search for collective nucleon flow in 
heavy ion collisions (l5| . The close connection between 
flow and jets continues at RHIC, where we seek the rela- 
tion between minijets and "elliptic flow." 

Event-wise sphericity S is a measure of structure in 
single-particle density p(0, 4>) on the unit sphere. Wc 
use dyadic notation to reduce index complexity, analo- 
gous to vector notation Q m = ^2 i riu(ui(f>i). S (with 
r% — > Pi) describes a 3D quadrupole with arbitrary orien- 
tation. Given p = p [sin(0) cos(0), sin(0) sin(</>), cos(0)] = 
pu(9, 4>) we have 

n n 

2S = 2^2PiPi = 2^p 2 £(0 i , (C6) 

i=l i=l 
n 

= J2p 2 iU(6i,<t>i) 

i 

= n{p 2 U(e, (/>)), 
the last being an event-wise average, where 

U(e,<f>) =sm 2 (9)l + y(9,cj)) (C7) 

and 

y(M)= (C8) 

sin 2 (0)cos(2<?(>) sin 2 (0) sin(20) sin(20) cos(</>) " 
sin 2 (0)sin(2<?(>) - sin 2 (0) cos (20) sin(20)sin0 
sin(20) cos(0) sin(20)sin0 3cos 2 (0)-l 
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In terms of event-wise quadrupole moments q2 m derived 
from the Y^ m 



2S 



,(p 2 sin 2 (6»))J 



(C9) 




3q22 
-3tq 2 i 



3q22 
-3?q22 



-5Jq 2 i 



2 Q20 



an event-wise estimator of angular structure on the unit 
sphere, its reference defined by 2S re f — n(p 2 sin 2 (#)) X, 



and p 2 sin 2 (8) = 
normalized to S 
Note that 



pf. The sphericity tensor of [14( was 
= S/n(p 2 ). 



Q = 3S-n(p 2 )X 



(CIO) 



is the traceless Cartesian quadrupole tensor appearing in 
the Taylor expansion of the (p equivalent of the) electro- 
static potential [38] . We have defined instead 



0! = 35 



(Cll) 



an alternative quadrupole tensor wherein each element 
is a single spherical quadrupole moment. The difference 
lies in the diagonal elements: linear combinations of the 
3?q2 m in the diagonal elements of Q are simplified to 
single moments in Qf. The ensemble mean of both ten- 
sors for an uncorrelated (spherically symmetric) system 
or system with event-wise quadrupole orientations ran- 
domly varying is the null tensor (all elements zero). 



APPENDIX D: SUBEVENTS 

The "subevent" is a notional re-invention of partition- 
ing/binning, the latter having a history of more than a 
century in mathematics. In conventional flow analysis 
subevents are groups of particles in an event segregated 
on the basis of random selection, charge, strangeness, 
PID or a kinematic variable such aspt, y or rj. The scalar- 
product method [3(| is based on a covariance between two 
single-particle bins (nominally equal halves of an event). 
The subevent method is thus a restricted reinvention of 
a common concept in multiparticle correlation analysis: 
determining covariances among all pairs of single-particle 
bins at some arbitrary binning scale - a two-particle cor- 
relation function. Diagonal averages of such distributions 
are the elements of autocorrelations. 

In the language of conventional flow analysis one way 
to eliminate statistical reference Q 2 e j from Q 2 m is to par- 
tition events into a pair of disjoint (non-overlapping) 

subevents A, B [12]. In that case Q ma -Qmb = V ma -V mb = 
n anbV^ nab , a covariance. The partition may be asym- 
metric (unequal particle numbers) and may be as small 
as a pair of particles. In addition to eliminating the 
self-pair statistical reference such partitioning is said to 
reduce nonflow correlation sources, depending on their 



physical origins and the partition definition [30J. We as- 
sume for simplicity that there is no nonflow contribution. 
Subevent pairs can be used to determine the event-plane 
resolution for subevents A, B and full events. 

First, we consider the symmetric case, defining equiv- 
alent subevents A and B with multiplicities tia = tib = 
n/2 from an event with n particles. E.g., subevent A 
has azimuth vector Q m A = X)™gA u(m<t>i). The scalar 
product is a covariance 

Qma ■ Qmb = Qa Qh(cOs(m[* Q - $ 5 ])) (Dl) 

= ^2 cos(m[^j - <j>j]) 

= n a 71b Vmab = Knafe' 



with e.g. Q 2 a = n a + n a (n a - 1) 



cos(m[* ma - * m &]) 



r.2 



V ma- Then 



V 2 



=(D2) 



If subevents A and B are physically equivalent (e.g., a 
random partition of the total of n particles), then 



cos(m[* ma - fmb]) = r ab 



V 2 



V 2 



(D3) 



cos(m[^ ma - ^ r ]) cos(m[* mb - <Ev] 



where r ab = V 2 wb / VKL V mb is Pearson's normalized co- 
variance between subevents A and B for the m th power- 
spectrum elements. If A and B are perfectly correlated 
(r a b = 1) then 



cos(m[* ma - * m6 ]) = cos 2 (m[* ma - * r ]) (D4) 

In general, V^/n = (l + r<,&) Vma/^a, which provides the 
exact relation between the EP resolution for subevents 
and for composite events A + B. It is not generally cor- 
rect that cos(m[^f m — \f r ]) = \J~2 • cos(m[^ ma — ^r])- In 
this case 



cos 2 (m[*„ la - <J r ]) = 



V 2 



1 fla + V* 



(D5) 



and V^ a = for perfectly correlated subevents. 

Second, we consider the most asymmetric case A = 
one particle and B = n — 1 particles. 

(cos(m[<& - $ r ]) cos(m[* mi - * r ])) = (D6) 



y/n - 1 + (n - l)(n - 2)v, 



,~2 

III I 



n - 2 Q'm 



where Q„ = n — 1 + V^ 2 describes a subevent with n — 1 
particles. In general, the EP resolution for a full event of 
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n particles is given by 



cos 2 (m[* m - $V]) 



mi 



(D7) 



Measurement of the EP resolution is simply a measure- 
ment of the corresponding power-spectrum element, since 



COS 2 (m[W m - gV|) 

1 - cos 2 (m[* m - \E' r 



(D8) 



In [22J the approximation 

(cos(m[* m - vl/,,])) 2 « lv^/n (D9) 

is given for V^/n <C 1 or Q 2 „ ~ n. 

Equal subevents, as the largest possible event parti- 
tion, imply an expectation that only global (large-scale) 
variables are relevant to collision dynamics (e.g., to de- 
scribe thermalized events). The possibility of finer struc- 
ture in momentum space is overlooked, whereas autocor- 
relation studies with finer binnings and the covariances 
among those bins discover detailed event structure highly 
relevant to collision dynamics. 



and the most peripheral 10 or 20% bins therefore have 
large systematic errors resulting in abandonment. Flow 
measurements with EP estimation are also excluded from 
peripheral collisions due to low event multiplicities. 

In contrast, with the power-law method running inte- 
grals of the Glauber parameters and n c h can be brought 
into asymptotic coincidence for peripheral collisions re- 
gardless of the uncertainty in the total cross section. Pa- 
rameter v measures the centrality and greatly reduces 
the cross-section error contribution. Centrality accuracy 
< 2% on v is thereby achievable down to N-N collisions. 
That capability is essential to determine the correspon- 
dence of A-A quadrupole structure in elementary colli- 
sions, to test the LDL hypothesis for instance: is there 
"collective" behavior in N-N collisions? 

For central collisions the upper half-maximum point 
on the power-law minimum-bias distribution provides a 
precise determination of b = on n c h and therefore v. 
The 6 = point is critical for evaluation of correlation 
measures relative to Glauber parameter e in the context 
of hydro expectations for v%/ 'e. 



1. Geometry parameters and azimuth structure 



APPENDIX E: CENTRALITY ISSUES 

Accurate A-A centrality determination and the cen- 
trality dependence of azimuth multipoles and related pa- 
rameters is critical to understanding heavy ion collisions. 
We must locate b = accurately in terms of measured 
quantities to test theory expectations relative to hydro- 
dynamics and thermalization. And we must obtain accu- 
rate measurements for peripheral A-A collisions to pro- 
vide a solid connection to elementary collisions. 

1. Centrality measures 

In [33] is described the power-law method of centrality 
determination. Because the minimum-bias distribution 
on participant-pair number n part /2 goes almost exactly 
as (n part /2)~ 3 / 4 the distribution on (npart/2) 1 / 4 is al- 
most exactly uniform, as is the experimental distribution 
on rij^ dominated by participant scaling. Those sim- 
ple forms can greatly improve the accuracy of central- 
ity determination, especially for peripheral and central 
collisions. The cited paper gives simple expressions for 
Tipart/2, ribi n and v relative to fraction of total cross sec- 
tion. 

In conventional centrality determination the minimum- 
bias distribution on n c h is divided into several bins rep- 
resenting estimated fractions of the total cross section. 
The main source of systematic error is uncertainty in 
the fraction of total cross section which passes triggering 
and event reconstruction. The total efficiency is typi- 
cally 95%, the loss being mainly in the peripheral region, 



We consider the several A-A geometry parameters rel- 
evant to azimuth structure. In Fig. [TT] (left panel) we 
plot rip ar t/2 vs v using the parameterization in [32| . 
The relation is very nonlinear. The dashed curve is 
n par t/2 ~ 2^ 2 57 . The most peripheral quarter of the 
centrality range is compressed into a small interval on 
n P art/2. Mean path-length v is the natural geometry 
measure for sensitive tests of departure from linear N- 
N superposition, whereas important minijet correlations 
(nonflow) oc v are severely distorted on n par t- 
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FIG. 11: Left panel: Participant pair number vs mean path- 
length v for 200 GeV Au-Au collisions. Because of the nonlin- 
ear relation the peripheral third of collisions is compressed to 
a small interval on n vart /2. Right panel: Impact parameter b 
vs v. To good approximation the relation is linear over most 
of the centrality range. 

In Fig. [Tl] (right panel) we plot impact parameter b vs 
v, again using the parameterization in [32j with fractional 
cross section <t/<tq = (b/b ) 2 . We note the interesting 
fact that over most of the centrality range b/bo ~ (R — 
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v)/(R - 1), with b = 2R= 14.7 fm for Au-Au. Thus, 
any anticipated trends on b are also accessible on v with 
minimal distortion. 

In Fig. [T2l (left panel) we show the LDL parameter 
1/ S dn c h/dri [2J] vs v for three collision energies. The 
energy dependence derives from the multiplicity factor, 
which we parameterize in terms of a two-component 
model [39( | . Weighted cross-section area S(b/b ) (fm 2 ) 
is an optical Glauber parameterization from [3l| . Both v 
and 1/ S dn^/dr] are pathlcngth measures. They can be 
compared with the inverse Knudson number K~ x intro- 
duced in [3|| as a measure of collision number. The LDL 
measure is based on energy-dependent physical particle 
collisions, whereas v is based on A-A geometry alone. 
The relation is monotonic and almost linear. Thus, struc- 
ture on one parameter should appear only slightly dis- 
torted on the other. 




10 20 
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FIG. 12: Left panel: Correspondence between LDL parame- 
ter 1/S dn c h/drj and centrality measure v for three energies. 
Right panel: Theory expectations for two limiting cases at 
200 GeV. The solid curve is derived from the solid curve in 
Fig. 1101 (upper-right panel) using the relation in the left panel. 

The hatched region is typically not measured in a con- 
ventional flow analysis, due to a combination of large sys- 
tematic uncertainty in the centrality determination and 
large biases in flow measurements due to small multi- 
plicities. However, peripheral collisions provide critical 
tests of flow models: e.g., how does collective behavior 
(if present) emerge with increasing centrality? In this 
paper we describe analysis methods which, when com- 
bined with the centrality methods of [32], make all A-A 
collisions accessible for accurate measurements down to 
N-N 



In Fig. [12] (right panel) we show t>2/e vs 1 / S dn c h / drj 
for theory expectations (hatched bands) and the simula- 
tion in Sec. IVII CI The latter is based on a simple error 
function on v and is roughly consistent with four-particle 
cumulant results at 200 GeV [2||. We observe that the 
solid curve is not consistent with either the LDL trend 
for peripheral collisions (the LDL slope is arbitrary) or 
the hydro trend for central collisions. That provocative 
result suggests that accurate analysis of azimuth corre- 
lations over a broad range of energies and centralitics 
with the methods introduced in this paper and [32| may 
produce interesting and unanticipated results. 

3. Correlation measures 



If the centrality dependence of azimuth structure is 
to be accurately determined the correlation measure em- 
ployed must have little or no multiplicity bias, includ- 
ing statistical biases and irrelevant multiplicity factors 
which lead to incorrect physical inferences. The quantity 
Ap/y/p re f is the unique solution to a measurement prob- 
lem subject to multiple constraints. It is the only portable 
measure (density ratio) of two-particle correlations ap- 
plicable to collision systems with arbitrary multiplicity. 
Ap/^/p re f is invariant under linear superposition. If, ac- 
cording to that measure, central Au-Au is different from 
N-N the difference certainly indicates a unique physical 
aspect of Au-Au collisions relative to N-N, exactly what 
we require in a correlation measure. Conventional flow 
measures do not satisfy that basic requirement. 

Drawing a parallel with measures of (pt) fluctuations 
we compare V2 



-'pt 



. Both are square roots of 
per-pair correlation measures which tend to yield mis- 
leading systematic trends (on centrality and energy) [Io| . 
In contrast V£ «-> AE 2 t .„ [5j (the total variance dif- 
ference for p t fluctuations) are integrals of two-particle 
correlations (azimuth number correlations vs pt correla- 
tions). The first is a measure of total azimuth corre- 
lations, the second a measure of total pt variance dif- 
ference [4l[, the integral of a two-particle distribution 
relative to its reference. In a minimally-biased con- 



text we then have v n 



V^ n /n{n — 1), analogous to 



S Pt (CERES) = ^AE% t .Jn(n- l)p 2 t [3 (the similar- 
ity of notations is coincidental) as the square roots of 
per-pair measures. Both are physically misleading. 
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